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SECTION  I.  INTRODUCTION 
A.  Study  Objective 

The  purpose  of  this  study  is  to  perform  research  necessary  to 
develop  new  discrete  design  techniques,  discrete  adaptations  of  control 
theory,  and  in  particular  extended  Parameter  Space  analysis  techniques 
that  will  help  to  insure  development  of  effective  guidance  and  control 
(G&C)  systems  for  future  weapon  systems.  This  research  will  permit  the 
development  of  candidate  advanced  G&C  systems  for  the  family  of  Future 
Army  Modular  Missiles  (FAMM).  The  prime  output  of  this  technical  effort 
is  to  be  the  enhancement  of  the  ability  to  create  and  evaluate  new  discrete 
G&C  systems  capable  of  meeting  tht  elusive  threat  of  the  future. 

The  study  effort  is  composed  of  seven  tasks  described  in  the  Scope 

of  Work  (as  altered  by  the  three  modifications  to  the  original  contract), 

contained  in  Technical  Requirement  No.  T-0208  entitled  "Application  of 
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Digital  Technology  to  Guidance  and  Control  Theory."  The  requirements 
imposed  by  these  seven  tasks  are  contained  in  the  referenced  contract  docu¬ 
ments  and  re-stated  (for  convenience  of  the  reader)  below. 

Task  I  requirements: 

"2.1  Investigate  to  determine  if  a  sensitivity  problem  due  to 
clustering  of  the  system  poles  and  zeros  exists  for  the  parameter  plane 
technique  of  digital  system  analysis/design." 

"2.2  If  a  sensitivity  problem  (i.e.,  clustering  of  the  system 
poles  and  zeros)  is  found  in  the  investigation  of  Paragraph  2.1,  a  technique 
shall  be  developed  to  mathematically  spread  apart  the  system  poles  and  zeros 
via  intermediate  complex  transform  techniques  to  eliminate  the  sensitivity." 

Task  II  requirements  : 

"2.3  Utilizing  Government  supplied  data  (Section  IV,  page  19,  of 
contract,  reference  1)  on  such  subjects  as,  but  not  limited  to,  projected 
future  targets,  missile  sensor-tracker  technology,  airframe  performance 
capabilities/limitations,  onboard  microprocessor  technology,  and  predicted 
computer  analysis  techniques,  the  contractor  shall  propose  at  least  two 
candidate  Advanced  Theoretical  Optimal  Guidance  and  Control  Systems." 
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Task  III  requi rements  : 

"2.7  The  contractor  shall  investigate  the  capabilities  and 
limitations  of  utilizing  a  relatively  new  method  of  analysis  and  design  of 
digital  systems  called  the  Parameter  Space  technique  to  analyze  guidance 
and  control  systems  and  shall  extend  and  refine  the  Parameter  Space  technique 
through  efforts  to  reduce  the  number  and  severity  of  the  limitations  so  noted. 
Also,  other  digital  design  techniques  shall  be  expanded  and  restructured  to 
make  them  usable  to  the  applications  engineer.  This  shall  include,  but  not 
be  limited  to,  investigating  the  possibility  of  not  having  to  obtain  the 
system  characteristic  equation  before  applying  the  Parameter  Space  technique." 

"2.8  The  digital  system  design  and  analysis  techniques  considered 
in  sub-section  2.7  shall  be  fully  documented  with  at  least  one  example  for 
aid  in  applying  the  techniques." 

Task  IV  requirements : 

Deleted  (per  Mod.  P00002).3 

Task  V  requirements: 

"3.1  The  contractor  shall  conduct  a  technical  seminar  for  selected 
members  of  the  Guidance  and  Control  Directorate.  The  purpose  of  the  seminar 
is  to  develop  a  design  expertise  among  members  of  the  Directorate  to  enable 
them  to  design  and  evaluate  digital  (discrete)  guidance  and  control  systems. 
Both  existing  and  newly-developed  sampled-data  techniques  will  be  described 
and  discussed  in  detail  with  practical  applications  conducted  during  the 
course  of  the  seminar  conferences." 

"3.2  The  contractor  shall  document  a  comprehensive  outline  of 
the  material  presented  in  the  seminar  defined  in  sub-section  3.3  as  specified 
on  the  control  data  documents  list  (DD  Form  1423)." 

Task  VI  requirements: 

"3.3  The  contractor  shall  apply  unique  digital  design  and  analysis 
expertise  and  experience  to  evaluate  planned  and  proposed  digital  guidance  and 
control  systems  for  the  PERSHING  II.  When  it  is  deemed  appropriate  (subject 
to  Government  approval),  the  contractor  will  recommend  modifications  to  these 
systems  with  the  aim  of  obtaining  Improved  system  performance  and/or  de¬ 
creased  system  complexity/cost.  To  analyze  promising  systems  the  contractor 


will  (subject  to  Government  approval)  either  conduct  simulations  or  recom¬ 
mend  their  conduct  by  the  Guidance  and  Control  Analysis  Group  (DRSMI-RGN) 
"3.4  The  contractor  shall  document  the  analyses  and  evaluations 
indicated  in  sub-section  3.3  as  specified  on  the  control  data  requirement 
list  (DO  Form  1423)." 

Task  VII  (Unnumbered  in  contract)  requirements: 

"a.  Investigate  present  digital  filter  requirements  and  recast 
them  for  the  desired  performance  rather  than  emulating  analogue  forms." 

"b.  Identify  technical  areas  to  be  expanded  or  initiated  in 
order  to  develop  future  digital  filters  more  effectively  and  efficiently." 

B.  Study  Schedule 

The  requirements  of  Task  I  were  completed  on  30  September  1979. 

They  were  documented  in  a  technical  report  entitled,  "A  Mathematical  Means 
of  Spreading  Clustered  System  Roots. The  requirements  of  Tasks  1 1- VI I 
are  to  be  completed  by  30  September  1980. 

C.  Contents  of  Report 

This  technical  report  covers  the  efforts  performed  in  completing 
the  requirements  of  Tasks  1 1- VI I  of  the  study.  As  stated  in  the  Documenta¬ 
tion  Addendum  to  Form  1423  as  amended  by  contract  modifications,  the  con¬ 
tractor  is  required  to  document  the  following  tasks:  Task  I,  Task  III, 

Task  V,  and  Task  VI. As  indicated  in  Part  B  (above)  of  Section  I,  the 
documentation  requirement  for  Task  I  has  been  met  by  the  submission  of  a 

5 

technical  report.  Although  not  required,  a  number  of  efforts  that  were 
performed  under  Task  II  are  documented  herein.  Also,  although  not  required, 
a  summary  of  the  work  performed  under  Task  VII  is  included. 
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SECTION  II.  ADVANCED  G&C  SYSTEMS  (Task  II) 


A.  Discussion 

The  discharge  of  this  task  has  been  centered  largely  around  the 
efforts  to  develop  and  implement  a  computer  simulation  that  is  sufficiently 
flexible  to  assess  and  compare  performance  of  various  candidate  G&C  systems. 
Additionally,  effort  has  been  expended  in  determining  the  G&C  system  aspects 
associated  with  the  Anti-Tactical  Ballistic  Missile  (ATBM)  field.  Finally, 
certain  efforts  were  expended  to  better  define  advanced  G&C  system  require¬ 
ments  and  potential  characteristics. 

B.  Computer  Simulation  Development 

The  development  of  a  suitable  computer  simulation  of  a  FAMM  system 
has  been  the  subject  of  numerous  meetings  with  personnel  of  the  G&C  Analysis 
Group  and  from  the  Aeroball istics  and  Propulsion  Directorates.  An  unclassified 
design  point  for  the  missile  at  the  terminus  of  its  trajectory  has  been  se¬ 
lected  as  Mach  8,  70,000  feet  altitude,  and  lateral  acceleration  of  12  g's. 

This  is  to  meet  a  threat  that  is  aerial  in  nature  and  combines  high  speed 
and  high  agility  (i.e.,  high  g-maneuver  capability)  with  low  radar  cross  sec¬ 
tion.  It  is  this  threat  that  leads  to  the  requirement  for  developing  an 
optimally  guided  and  controlled,  highly  maneuvering,  defensive  missile  using 
terminal  guidance  sensors  chosen  from  across  a  wide  range  of  the  frequency 
spectrum.  The  attendant  required  high  level  of  control  authority  demands 
that  the  airframe  and  propulsion  be  developed  beyond  the  present  state  of  the 
art:  hence  the  involvement  of  Aeroball istics  and  Propulsion  Directorate 
personnel  in  the  above  referenced  discussions.  However,  in  order  to  have  a 
simulation  that  is  checked  out  and  running,  a  simplified  model  (termed  FAMM  I) 
has  been  arrived  at.  It  will  initially  model  the  airframe  as  a  standard  one 
that  already  exists  and  is  coded.  The  actuators  and  sensors  will  be  added  as 
program  modules  and  will  initially  have  dynamic  transfer  functions  of  unity. 

An  existing  Proportional  Navigation  Guidance  (PNG)  module  will  be  included 
with  the  option  of  later  replacing  it  with  an  existing  optimal  G&C  module. 

The  foregoing  systems  are  well  documented  for  the  G&C  Analysis  Group  by  M&S 
Computing  Company,  so  they  are  not  repeated  in  this  report. 
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In  August  1980,  it  was  decided  to  begin  investigating  the  appli¬ 
cability  of  the  advanced  G&C  systems  being  studied  to  the  demanding  pro¬ 
blems  posed  by  the  new  Anti-Tactical  Ballistic  Missile  (ATBM)  systems. 

The  threat  posed  by  enemy  TBM's  is  significantly  more  severe,  particularly 
from  a  dynamical  point  of  view,  than  that  used  in  developing  the  design 
point  of  Section  I  IB  (above). 

As  a  first  step,  on  8  July  1980,  Dr.  Pastrick  and  Dr.  Seltzer 
were  exposed  to  a  review  of  ATBM  G&C  that  is  under  investigation  by  the 
US  Army  Ballistic  Missile  Defense  Systems  Command  (BMSCOM).  This  was  follow¬ 
ed  by  a  visit  to  The  Charles  Stark  Draper  Laboratory,  Cambridge,  Massachu- 
ssetts ,  on  24  July  by  Dr.  Pastrick,  Mr.  Roy  Pugh  (G&C  Directorate),  and 
Dr.  Seltzer.  The  work  being  done  on  ballistic  missile  defense  systems  for 
BMSCOM  was  presented  in  detail  (Appendix  A). 

Shortly  after  effort  was  initiated  to  study  this  aspect  of  the 
problem.  Control  Dynamics  Company  support  was  provided  in  a  review  of  the 
US  Navy's  High  Angle  Threat  (HAT)  program.  A  review  of  that  program  was  con 
ducted  and  the  unclassified  portion  provided  Dr.  Pastrick  (Appendix  B).  As 
a  result  of  that  effort,  the  Control  Dynamics  Company  was  requested  to  con¬ 
tinue  participation  in  the  continuing  evaluation  of  the  HAT  program  with  the 
selection  of  the  most  promising  concepts  to  carry  into  any  further  studies. 
The  critical  issues  associated  with  the  HAT  defense  have  been  identified; 
they  are  shown  in  Appendix  C. 

A  basis  has  been  developed  from  which  a  future  study  of  the  ATBM 
implications  on  advanced  G&C  may  be  studied.  Further,  the  characteristics 
of  such  a  G&C  system  may  now  be  developed,  should  such  an  effort  be  desired. 

D.  Documentation 

There  are  some  items  of  documentation  that  are  not  referred  to  in 
other  portions  of  this  report.  They  are  described  in  this  section. 


1 .  Advanced  Analysis  for  Future  Missiles  (MI COM  TR) 

Dr.  Seltzer  organized,  reviewed,  and  proof-read  US  Army  Missile 
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Command  Technical  Report  RG-80-8,  "Advanced  Analysis  for  Future  Missiles," 
dated  21  November  1979.®  It  describes  the  work  that  was  completed  on  the 
(then)  new  program  to  develop  an  advanced  G&C  system  for  Future  Army 
Modular  Missiles.  It  includes  the  works  of  members  of  the  G&C  Analysis 
Group,  the  Aeroballistics  Directorate,  and  several  technical  specialists 
on  contract  from  universities  and  private  industry. 

2.  Future  US  Army  Missile  G&C  Systems  (AIAA  Paper) 

Dr.  Seltzer  prepared  and  presented  an  invited  (by  the  American 
Institute  of  Aeronautics  and  Astronautics)  paper  entitled,  "Future  US  Army 
Missile  G&C  Systems."  It  was  presented  at  the  1980  AIAA  International 
Conference,  Baltimore,  MO,  on  5-8  May  1980  (Appendix  D).^ 

3.  Special  Edition  of  AIAA  Journal  of  Guidance  and  Control 

Dr.  Pastrick  and  Dr.  Seltzer  were  invited  by  the  Editor-in- 
Chief  of  the  AIAA  Journal  of  Guidance  and  Control  to  organize  a  special 
edition  of  that  journal.  The  edition  will  be  devoted  primarily  to  tactical 
missile  G&C.  They  are  also  preparing,  with  Professor  M.  Warren  of  University 
of  Florida,  a  survey  paper  for  that  edition. 


There  are  some  trips  and  meetings  that  are  not  referred  to  in  other 
portions  of  this  report.  They  are  described  in  this  section. 

1 .  Guidance  Test  Symposium 

On  10-12  October  1980,  Dr.  Seltzer  attended  and  participated  in 
the  Guidance  Test  Symposium  at  Holloman  Air  Force  Base,  New  Mexico. 

2.  Office  of  the  Under  Secretary  of  Defense  for  Research  and 

Engineering 

On  19-20  November  1979,  Dr.  Seltzer  traveled  to  Washington,  D.C. 
to  meet  with  Dr.  Pastrick  (then  on  one  year's  duty  with  the  Office  of  the 
Under  Secretary  of  Defense  for  Research  and  Engineering)  on  missile  G&C 
matters.  While  there  he  also  participated  in  a  meeting  with  Mr.  Charles 
Bernard,  the  Director  of  Land  Warfare,  DDRE,  to  discuss  future  tactical 
missile  G&C. 

3.  Joint  Technology  Common  Intercept  Missile 

On  15  January  1980,  Dr.  Seltzer  attended  a  briefing  and  discussion 
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of  a  potential  Joint  Technology  Common  Intercept  Missile  held  at  the 
Advanced  Systems  Concepts  Office  (ASCO).  The  main  briefings  were  presented 
by  Mssrs.  Ben  Glatt  (NAVAIR),  Tom  Hamilton  (NAVWEPS),  and  Charles  Weaver 
(NAVSEA).  The  ASCO  personnel  responded  with  briefings  on  their  new  Army 
concept  for  an  Evolutionary  Surface- to-Air  Missile  system.  The  briefings 
were  followed  by  a  general  discussion  of  system  requirements  and  delinea¬ 
tion  of  possible  areas  of  cooperation  or  joint  development  efforts. 

4.  AFTL/DLO,  Eglin  Air  Force  Base 

On  29  April  1980,  Dr.  Seltzer  visited  the  8th  AFATL/DLO,  Eglin 
Air  Force  Base,  Florida  with  Dr.  Pastrick.  Discussions  took  place  and 
design  and  program  information  concerning  progress  and  plans  for  the  US 
MICOM  and  the  USAF  advanced  G&C  systems  were  exchanged.  Key  USA F  partici¬ 
pants  were  Captain  Riggs  and  Lieutenants  McLemoon  and  Vergez. 

5.  1980  Tactical  Missiles  Conference 

On  30  April  -  1  May,  Dr.  Pastrick  and  Dr.  Seltzer  attended  the 
1980  Tactical  Missiles  Conference  at  Eglin  Air  Force  Base. 
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SECTION  III.  PARAMETER  SPACE  TECHNIQUE  (Task  III) 

A.  Accompl ishments 

The  Parameter  Space  method  is  a  technique  for  determining  stability 
and  dynamic  characteristics  of  a  control  system  in  terms  of  several  selected 
system  parameters.  The  development  and  use  of  the  technique  Is  described  in 
US  Army  Missile  Command  Technical  Report  T-79-64,  entitled,  "Sampled-Data 

O 

Analysis  in  Parameter  Space,"  dated  June  1979  (Appendix  E).  In  essence, 
the  method  provides  an  analytical  tool  for  use  in  control  system  design  and 
analysis.  Although  not  necessary,  its  application  is  facilitated  by  augment¬ 
ing  the  analytical  results  with  graphical  portrayals  in  a  selected  multi- 
parameter  space.  In  its  present  form,  the  method  requires  that  the  controlled 
system  be  described  by  a  characteristic  equation  which,  for  sampled-data  or 
digital  systems,  may  be  expressed  in  the  z-domain.  The  method  enables  the 
designer  to  evaluate  graphically  (or  analytically)  the  locations  of  roots  of 
the  characteristic  polynomial.  Hence,  he  may  design  the  control  system  in 
terms  of  his  selected  performance  criteria  (also  known  as  "pole  placement"). 

He  is  able  to  see  the  effect  on  the  roots  of  changing  numerical  values  of 
several  system  parameters. 

Also  developed  under  an  earlier  contract  with  MICOM  is  the  Cross- 

Multiplication  method  for  determining  a  digital  control  system's  response  to 

a  deterministic  input.  This  method  is  described  In  US  Army  Missile  Research 

Command  Technical  Report  T-79-58,  entitled,  "Determination  of  Digital  Control 

g 

System  Response  by  Cross-Multiplication,"  dated  29  May  1979.  This  technique 
assumes  that  the  closed-loop  transfer  function  is  available  in  the  z-  or 
modified  z-transform  domain.  The  Real  Translation  Theorem  is  then  applied  to 
the  result,  yielding  a  difference  equation  in  the  time-domain.  This  may  be 
solved  for  the  system  response  in  terms  of  the  reference  (or  other)  input(s) 
to  the  system  as  well  as  in  terms  of  system  state  initial  conditions.  The 
technique  also  permits  the  determination  of  intra-sampling  responses  of  the 
system  by  using  either  the  Submultiple  Method  or  the  Modified  Z-Transform 
Method.  An  additional  strength  of  the  method  is  that  the  input  need  be  spec¬ 
ified  only  at  the  sampling  instants.  Hence  it  is  not  necessary  that  the 
input  be  described  in  a  continuous  (with  respect  to  time)  manner,  nor  is  it 
necessary  to  obtain  the  z-transform  of  the  input(s). 


8 


1 


The  Parameter  Space  method  and  Cross-Multiplication  technique 
have  each  been  programmed  for  usage  in  the  Hewlett-Packard  Series  9835A 
digital  computer.  The  program  language  is  an  advanced  version  of  BASIC. 

A  listing  of  the  programs  is  included  in  Appendices  F  and  G.  To  use  the 
Parameter  Space  program,  one  must  have  the  characteristic  equation  in  the 
z-domain.  The  program  then  provides  contours  that  portray  the  stability 
boundary  in  terms  of  a  selected  two-parameter  space.  It  then  determines 
and  plots  contours  corresponding  to  the  complex  and  real  roots  of  the 
system  characteristic  equation,  thereby  permitting  design  by  pole  placement. 
The  program  contains  the  option  of  being  able  to  determine  the  system  roots 
at  any  location  in  the  selected  parameter  space  (see  subroutine  entitled 
"SILJAK"  for  an  efficient  root  determination  routine). 

The  Cross-Multiplication  program  requires  that  the  desired  closed- 
loop  transfer  function  be  available  in  the  z-domain.  Of  course,  this  can 
be  determined  from  the  characteristic  equation  used  for  the  Parameter  Space 
program.  From  this  transfer  function,  the  output  can  be  provided,  at  the 
sampling  instants,  in  terms  of  any  desired  input  known  (deterministically) 
at  the  sampling  instants.  Several  typical  inputs  (unit  step,  ramp, 
sinusoidal)  are  provided  by  the  program  itself. 

The  Parameter  Space  Program  has  been  structured  so  that,  at  a 
later  date,  stability  regions  and  maps  of  system  roots  can  be  determined 
and  plotted  for  continuous-data  systems  in  the  Laplace(s)  domain.  This  option 
will  be  useful  to  provide  a  reference  or  comparison.  Also,  often  contemporary 
control  system  designers  (sadly)  perform  their  design  analysis  in  the  s-domain 
and  then  "merely  sample  the  system  fast  enough."  This  modification  will  per¬ 
mit  an  analysis  of  such  systems. 

Because  many  designers  prefer  to  design  in  the  w-domain  (bilinear 
transformation  from  the  z-domain),  the  Parameter  Space  Program  also  has  been 
structured  to  permit  its  inclusion  at  a  future  date.  This  is  contingent  on 
future  investigations  of  the  applicability  of  the  Parameter  Space  method  to 
w-domain  analysis. 
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B.  Example 

1 .  System  Description 

To  better  demonstrate  the  capabilities,  limitations,  and  general 
use  of  these  two  techniques  and  their  implementation  into  two  digital  pro¬ 
grams,  an  example  is  provided.  For  pedagogical  purposes,  it  represents  a 
simplified  aerospace  vehicle  under  control  of  a  digital  computer.  See 
Figure  1  (Section  X)  for  a  description  of  the  system  in  block  diagram  form. 
Capital  letters  are  used  to  denote  matrices,  and  capital  letters  with  under¬ 
bars  denote  vectors. 

2.  Application  of  SAM 

Suppose  the  design  goal  is  to  specify  a  digital  control  system 
that  will  result  in  a  stable  system  whose  performance  is  desirable  (in  some 
sense).  The  approach  to  be  used  will  be  to  obtain  a  mathematical  description 
of  the  system  dynamics  that  can  subsequently  be  simplified.  The  intermediate 
objective  is  to  obtain  a  closed- loop  transfer  function  in  the  z-domain.  The 
ultimate  objective  is  to  use  this  closed-loop  transfer  function  to  assess 
system  performance,  including  stability. 

Referring  to  Figure  1,  one  is  able  to  write  the  "System  Equa¬ 
tions"  (first  column  in  Table  I,  below)  describing  the  system  dynamics. 
Applying  the  "SAM"  (Simplified  Analytical  Method)  technique,  one  may  then 
obtain  the  "Modified  System  Equations"  (Column  2)  and  the  "Pulsed  Equations 
(Column  3)." 


Table  I.  SAM 


System  Equations 

Modified  Sys.  Eqns. 

Pulsed  Equations 

* =  GpFGhoT 

1  -  GpFGh<>“* 

i*  ■  eyws* 

(1) 

U  =  DE  * 

U  =  DE* 

U*  .  d*e* 

(2) 

E  =  R  -  Y 

E  =  R  -  Y 

E*  =  R*  -  Y* 

(3) 

Y  *  CX 

I "  CGp' V 

Y*  -  (CGpFGh0)*U* 

(4) 
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Now  any  desired  closed-loop  transfer  function  may  readily  be  obtained  from 
the  Table  I  equations.  For  the  purposes  of  this  example,  suppose  one  wishes 
to  determine  the  response  of  the  system  output  vector,  X,  in  terms  of  the 
system  input  vector,  R.  The  value  of  these  two  vectors  at  the  sampling  in¬ 
stants  are  denoted  by  the  star  operator  (*),  i.e.,  X*  and  _R*.  From  the 
third  ("Pulsed  Equations")  column,  one  readily  obtains  that  relation  with 
minimal  mathematical  manipulation: 


E*  a  R*  -  Y* 


=  **  -  <CGpFShO 


)*U* 


a  R*  -  (CGpFGh0)*D*E* 

«  [I  +  (CGpFGh0)*D*]_1  R*  . 


(5) 


This  expression  for  £*,  in  terms  of  R*,  may  now  be  substituted  into  Equation 
(2)  to  obtain  the  control  vector  U*  in  terms  of  R*: 


U*  a  D*E* 


■  D*  [I  +  (CGpFGh0)*D*]_1  R*  . 


(6) 


Finally,  this  expression  for  JJ*  may  be  substituted  into  Equation  (1)  to 
obtain  the  sought-after  expression  for  X*  in  terms  of  R*: 


r  - 

*  (GpFGh0)*D*  Cl  +  (CGpFGh0)*D*J_1  R* 


(7) 


In  general,  pulsed  quantities  may  be  expressed  directly  in  the  z-domain, 
i.e., 


3 


{X*}  -  X(z) 


(8) 
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Taking  the  z-transfcrm  of  each  side  of  Equation  (7)  leads  to 


X(z)  "  GpFGhO(z)  0{z)  CI  +  CGpFGhO(z)  0(z)r1  ^(z)‘  (9) 

where  the  following  notation  has  been  used: 


} 


* 

t 


\ 


\ 


3(VW*;  *'  yW2> 


(9a) 


and 


3 «  VW*1  CGPFGho(z)  •  <»> 

As  previously  stated,  the  expression  of  Equation  (9)  will  be  simplified  for 
pedagogical  purposes.  It  is  therefore  assumed  that  only  planar  rotational 
dynamics  of  a  rigid  body  will  be  investigated.  If  the  system  output  X  is 
assumed  to  be  the  time  rate  of  change  of  attitude,  i.e.,  9,  then  the  actuator 
output  V  (assumed  scalar  in  form)  is 

V  *  19,  (10) 

where  I  represents  the  principal  moment  of  inertia  along  the  axis  of  rotation. 

In  the  Laplace  domain,  this  leads  to  the  definition  of  the  plant  transfer 
function,  Gp(s) 


If  the  zero  order  hold  is  assumed,  then 


(11) 


I 


t 

i 


Ws) 


and 


(12) 


12 


•ffVtririX 


VW2)  *3{Vs)  Ws)} 


*1 


1  -  e' 

"77 


•sT 


nbn 


03) 


Suppose  one  assumes,  for  the  form  of  the  digital  controller,  that 
olde  standby,  PID.  In  analogue  form,  the  controller  would  be 

D(s)  *  ^  +  ^  +  <D  • 

If  one  applies  the  trapezoidal  approximation  for  integration,  the  s”^  is 
replaced  by 


(14) 


T 

2 


This  substitution  into  Equation  (14)  leads  to  the  z-form  of  the  controller, 
i  .e. , 


D(z) 


+ 


(z+1)2 

(z-1)2 


+  Kr 


(15) 


For  simplification,  let 


and 


d.  p 
'  21  ’ 

(16P) 

K  T2 

d.  V 

*  41  ’ 

(161) 

j.  V. 

( 16D) 
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Equations  (15)  -  (16)  lead  to  the  expression, 


TD (z)  .  («p  ♦  ai  +  aQ)  z  +  2  (ar  -  ap)  z  +  (-ap  +  ar  +  ap) 
1 

From  Equation  (13),  one  may  obtain  the  expression 


IGo<W2>  1 

T  "  z-1  • 


07) 


(18) 


Multiplying  Equations  (17)  and  (18)  together,  one  obtains  an  expression  for 
G  G^q(z)D(z)  for  use  in  Equation  (9),  i.e., 

SS0(«)"«)  ■  ^  V  Z2  *  2  ■  -  *  1  ^  («. 
p  hO  (2-1 y 

The  scalar  version  of  vector-matrix  Equation  (9)  leads  to  the  elementary 
closed  loop  transfer  function. 


(20)  j 

] 

9  q  I 

Letting  X(z)  =  e(z)  and  R( z )  =  ec(z)  and  assuming  that  F  =  1  and  C  s  1,  one  \ 

finally  obtains,  with  the  use  of  Equation  (19),  the  desired  form  of  the  j 

closed-loop  transfer  function:  j 


X(z)  GpF6h0(z)P(z) 

R(z)  '  1  +  CGpFGh0(z)D(z)  ’ 


e(z)  s  -T  h 

M*>  ’  iT  k  *  £  h  ’ 


(21) 


where 
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(22) 


A  »•  [(-ap  +  aj  *  aQ),  2  (aj  -  aQ) ,  (ap  +  aj  +  aD)]T  , 

£  d*  [-1,  3,  -3,  1JT  ,  (23) 


In  *  Cl-  z«  Z^T  » 

and 


(24) 


ZD  *•  [1,  Z,  Z2,  Z3]T  .  (25) 

Now  consider  Equation  (21).  It  is  in  a  form  to  apply  both  the  Parameter 
Space  and  the  Cross-Multiplication  techniques,  since  its  denominator  is, 
when  set  equal  to  zero,  the  system  characteristic  equation.  The  only  para¬ 
meter  that  needs  to  be  specified  numerically  are  ap,  a^,  and  aQ.  These  will 
determine  the  response  of  the  system  to  any  specified  input.  Notice  that  no 
knowledge  of  complex  variables  has  been  needed  by  the  system  designer:  only 
a  blind  reliance  on  some  simple  rules  (SAM).  At  this  point  one  may  (blindly) 
determine  the  system  response  (by  Cross-Multiplication)  without  applying  the 
Parameter  Plane  method.  However,  application  of  the  Parameter  Space  tech¬ 
nique  (or  some  other  stability  and  system  performance  method  or  methods)  will 
determine  system  stability  and  desired  performance  requirements.  If  done  in 
terms  of  the  three  parameters  ap,  a^,  and  a^,  one  will  determine  the  numerical 
values  for  digital  controller  gains  Kp,  Kj ,  (normalized  by  the  moment  of 
inertia  I)  and  the  sampling  period  T. 

3.  Application  of  Parameter  Space  Technique 

First  the  control  system  designer  must  specify  two  parameters 
of  interest.  Since  he  has  three  (ap,  aj,  a^)  at  his  disposal,  assume  he 
selects  ap  and  a^  (since  he  may  ultimately  choose  to  not  use  any  integral 
gain).  He  must,  however,  choose  a  numerical  value  (which  may  be  iteratively 
changed  in  the  design  sequel)  for  a^.  Let  us  assume  a  value  of  aj  *  0.2  for 
this  example.  Then  the  characteristic  equation  --  the  denominator  of  Equation 
(21)  --  may  be  re-cast  in  the  following  form: 
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(26) 


C.E.(z)  =  AZQ  PI  +<BZ[)  P2  + 

r*J 

F-zo 

=  0 

(26) 

where 

PI  i-  a„ 

(27) 

and 

P2  =•  ip  ■ 

(28) 

Grouping 

the  coefficients  of  the  characteristic  equation 

in  the 

form 

of 

Equation 

(26),  i.e., 

8TID 

+  4TZN  *  [-I,  3,  -3,  1]  [1,  Z, 

z2, 

Z3]T 

+  C(-ap  +  ar  +  aQ) ,  2  ( 

aI  - 

aD),  (ap  +  ax 

+  aD)]  [1, 

z,  z2f 

■  (1  -  2Z  +  Z2)  aD  +  (- 

1  + 

Z2)  ap 

+(-l  +  aj)  +(3+  2a^ ) 

z  + 

(-3  ♦  a,)  Z2 

+  Z3  . 

(29) 

By  association,  matching  Equation  (29) 

with 

Equation  (26) 

leads 

to 

7  = 

[A  (0),  A  (1),  A  (2),  A  (3)]  = 

[1,  ■ 

-2,  1,0], 

(30) 

7  = 

[B  (0),  B  (1),  B  (2),  B  (3)]  =» 

[-1, 

0,  1,  0]  , 

(31) 

and 

F  = 

CF  (0),  F  (1),  F  (2),  F  (3)3  = 

C(-l 

+  a^ ) ,  ( 3  +  2aj ) ,  (- 

-3  +  < 

»j),  1] 

s 

[-0.8,  3.4,  -2.8,  1]  . 

(32) 
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Recognizing 
is  ready  to 
exposition , 
development. 


that  the  order  of  the  characteristic  equation  is  three,  one 
use  the  parameter  space  program,  MICPP  (Appendix  F) .  For 
the  computer  printout  will  be  included  in  the  following 


Z  TRANSFORM  METHOD 
THE  SYSTEM  ORDER  IS  =  3 


ELEMENTS  OF  A  MATRIX  A  <  0  :• ,  A  <  l  >  ,  A  <  2  >  ,  .  .  . 

1 . 0000E+00  -2.O000E+00  1.0000E+00  0.0000E+O0 


ELEMENTS  OF  B  MATRIX  B<0.) ,  B<  1  )  ,  B<2)  ,  .  .  . 

"1.0  0  0  0  E  +  0  0  0.0000E  +  00  1.0000E+00  0.0000E+00 


ELEMENTS  OF  F  MATRIX  F  <  0  > ,  F  <  1  >  ,  F  <  2  ,  .  .  . 
-S.0000E-01  3.400OE+00  -2.S000E+00  1.000OE+00 


First  the  stability  boundaries  will  be  determined  (these  are 
usually  associated  with  the  unit  circle,  or  c  *  0  contour,  in  the  z-domain 
and  with  the  singularities  a  (z  =  +  1  and  z  *  -1).  The  parameter  plane 
boundary  associated  with  z  =  -1  is  found  to  be  PI  s  aQ  *  2,  and  the  boundary 
associated  with  z  =  +  1  is  P2  =  ap  =  0.  For  plotting  convenience,  the  x  and 
y  intercepts  of  these  straight  lines  are  also  indicated. 


ZETA  =  0.0000E+00 

THE  STABILITY  BOUNDARIES  FOR  Z*+l  AND  -1  ARE  AS  FOLLOWS 

FOR  Z»-l 

PI  «  2.0000E+00 

THE  AXES  INTERCEPTS  ARE  AS  FOLLOWS: 

AND  AT  P2*0  ,  PI"  2.0000E+00 
PI  «  CONSTANT  *  2.0000E+00 
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FOR  Z=+l 

COEFFICIENTS  OF  PI 


P2  EQUAL  TO  ZERO 


THE  AXES  INTERCEPTS  FIRE  FlS  FOLLOWS: 
PI  =  CONSTFlNT  =-9. 9999E+99 
P2  =  CONSTANT  =-9.9999E+99 


For  the  boundary  associated  with  ?  =  0,  it  is  specified  by  the  designer  to 
let  the  angle,  u)nT,  vary  between  zero  and  pi  in  20  equal  increments.  Print¬ 
ed  out  are  these  values  of  and  the  associated  values  of  Pi  and  P2.  To 
investigate  the  possibility  of  the  existence  of  the  fourth  stability  boundary 
associated  with  the  case  where  the  Jacobian  becomes  identically  equal  to 
zero,  the  Jacobian  is  printed  out  in  the  second  column.  If  it  were  to  change 
sign,  then  there  obviously  would  be  a  stability  boundary  at  the  value  of 
uinT  where  it  (the  Jacobian)  was  equal  to  zero.  Also,  the  sign  of  the  printed 
out  Jacobian  indicates  on  which  side  of  the  x,  -  0  contour  the  stability 
region  lies. 


INITIAL  OMEGA  =  0 . OQOOE+OQ  FINAL  OMEGA  =  3.141SE+00 

NUMBER  OF  STEPS*  100 
Wr.T  =  1.0000E-25 

THE  DETERMINANT  OF  THE  COEFFICIENT  MATRIX  IS  EQUAL  TO  ZERO 
MATRIX  AB  = 

0.0000E+00  0.0000E+O0 

0.0000E+00  2.0000E-25 


MATRIX  Fp  = 
-3. 0O00E-01 


-3. 00OOE-26 


WnT »  3. 1416E-02  DET--E .  1 997E-05 
WnT-  E.2332E-02  BET--4 . 95E 1 E-04 
WnT-  9.4243E-02  DET-- 1 . E706E-03 
WnT  =  1.25EEE-01  DET  =  -3 . 9332E-03 
WnT-  1 . 570SE-9 1  DET--7. 7039E-03 
WnT ■  1 . 8850E-0 1  DET  =  - 1 . 327EE-02 
WnT  =  2.1991E-01  BET--2 . 1 0 1 4E-02 
WnT-  2.5133E-91  DET  — 3.  1252E-62 
WnT  =  2.8274E-01  BET  —  4. 431  IE-02 
WnT-  3.1416E-01  DET--E.0497E-02 
WnT-  3.4538E-01  DET— 3 . 0 1 04E-02 
WnT-  3.7639E-01  DET--1 . 0340E-01 
WnT-  4.0841E-01  DET  =  - 1 . 3865E-0 1 
WnT-  4 . 3932E-0 1  DET— 1 . 6209E-0 1 
WnT  =  4.7124E-81  DET—  1 . 9793E-0 1 
WnT-  5.0265E-01  DET— 2. 3836E-0 1 
WnT-  5 . 3487E-0 1  DET  — 2. 335SE-0 1 
WnT-  5.  6549E-61  DET— 3. 3365E-81 
WnT-  5.9690E-01  DET  — 3. 3878E-01 
WnT-  S.2832E-01  DET— 4 . 4903E-01 
WnT-  E.  5973E-0 1  DET— 5.  1446E-01 
WnT-  S .  9 i  1 SE-0 1  DET  — 5.SS12E-01 
WnT-  7 . 22S7E-0 1  BET— E .  E 1 02E-0  1 
WnT-  7 . 539SE-0 1  DET  — 7. 4213E-31 
WnT-  7.3540E-01  DET--8. 2843E-01 
WnT-  8.1631E-01  DET— 9.  1982E-81 
WnT-  8 •  4823E-0 1  DET— 1 . 01 62E+03 
WnT-  3.79S3E-01  DET— 1 .  1 1 73E+00 
WnT-  9.11 08E-0 l  DET--1 . 2235E+00 
WnT  =  9 •  424SE-0 1  DET— 1 . 3340E+00 
WnT-  9 . 7389E-0 1  DET— 1 . 4433E  +  00 
WnT-  1 . 0053E+80  DET  —  1 . 5E77E+00 
WnT-  1.03E7E+00  DET-- 1 . 6904E+00 
WnT-  1.0631E+88  DET— 1 .8166E  +  00 
WnT  =  1 . 0996E+00  DET  «- 1 . 94E0E+00 
WnT-  1 . 1310E+00  DET--2. 0733E+00 
WnT-  1 . 1 E24E+8©  DET--2. 2131E+00 
WnT-  1 .  1 938E  +  00  DET— 2 . 3500E+00 
WnT-  1 . 2252E  +  00  DET— 2. 4387E+00 
WnT*  1 . 25EEE+00  DET--2 . 6287E+00 
WnT-  1.2881E+00  DET--2 . 7693E+00 
WnT-  1.31  93E  +  00  DET— 2 . 9 1  03E+00 
WnT-  1.3309E+80  DET  — 3 . 832 1E  +  00 
WnT-  1 . 3S23E  +  00  DET— 3.  1929E  +  0© 
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The  stability  boundaries  may  now  be  plotted  on  the  a^  -  ap  parameter  plane 
(Figure  2).  The  stable  region  lies  in  the  first  quadrant,  bounded  (on  the 
left)  by  the  c  ■  0  curve  and  (on  the  right)  by  the  PI  =  a^  *  2  line. 

Since  the  characteristic  equation  is  third  order,  it  must 
have  three  roots,  c-contours,  printed  out  as  functions  of  the  argument 
wnT,  display  the  characteristics  of  complex  conjugate  pairs  of  roots. 

Since  the  c  a  0  contour  exists,  there  must  be  some  complex  conjugate  roots. 
However,  since  there  are  three  roots,  there  also  must  exist  a  real  root. 
Hence,  at  each  point  in  the  a^  -  ap  parameter  plane,  one  must  expect  the 
existence  of  either  one  complex  conjugate  pair  of  roots  plus  one  real  root 
or  three  real  roots.  It  will  be  shown  that,  for  this  case,  the  former 
prevails. 

Zeta  contours  for  2eta  values  of  0.2,  0.5,  0.707,  and  i  are 
shown  as  functions  of  oinT  in  the  following  printouts..  These  contours  are 
also  plotted  on  Figure  2.  For  clarity,  the  associated  values  of  uinT  are 
omitted  from  the  figure. 
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1 . 8696E+00 

P2  = 

1.29426+00 

PI* 

1 . SS09E+0O 

P2  = 

1 . 30036+00 

PI* 

1 . 3902E  +  0O 

P2* 

1 . 30456  +  00 

PI* 

1 . 3976E+00 

P2= 

1 . 3O676+00 

PI* 

I . 9031 E+00 

P2  = 

1.30696+00 

PI* 

1 . 90S8E+0O 

P2* 

1.30536+00 

Pi* 

1 . 9OS7E  +  00 

P2  = 

1 . 30186+00 

Pl  = 

1 .9088E  +  00 

P2* 

l . 29656+00 

OMEGA  *  3 . 141 ££+0O 


22 


THE  DETERMINANT  OF  THE  COEFFICIENT  MATRIX  IS  EQUAL  TO  ZERO 
MATRIX  AB  ■ 

0.00006+00  0.0000E+00 

O.0000E+00  1.7321E-25 

MATRIX  Fp  =  .  . 

-8.0000E-01  -S .  9£s2E-*:o 


WnT  = 

6.23326-02 

BET*-4.0321E-04 

PI* 

2. 0264E+02 

P2* 

6. 3702E+00 

WnT* 

1 . 2366E-01 

DET*-3. 0232E-03 

PI* 

3. 0720E+O 1 

P2* 

3.  1 947E+O0 

WnT* 

1 . S83OE-0 1 

DET=-9. 5505E-03 

PI* 

2. 2637E+01 

P2* 

2.  1 445E  +  00 

WnT* 

2.51 33E-0 1 

DET  *-2 . 1 1626-02 

PI* 

1 . 2S49E  +  0  1 

P2* 

1 .6277E+00 

WnT* 

3.  1416E-01 

DET=-3 . S535E-02 

PI* 

3.351 3E  +  00 

P2  = 

1 . 3257E+00 

WnT* 

3. 7699E-01 

D6T*-6. 2161E-Q2 

PI* 

5 . 9353E+00 

P2* 

1 .  1 322E+O0 

WnT* 

4. 39326-01 

D  E  T  *  -  9 .  1 9  0  4  E  -  O  2 

PI* 

4 . 53306+00 

P2  = 

1 . 0O12E+00 

WnT* 

5.02656-01 

DET=-1 . 27566-01 

PI* 

3.3936E+00 

P2* 

9 . 0990E-0 1 

WnT* 

5. 63496-01 

DET =- 1 . 6364E-0 1 

Pi* 

2.9831E+30 

P2  = 

3 . 4333E-0 1 

WnT* 

6. 2832E-01 

DET—2.  1432E-01 

PI* 

2. 5703E+O0 

F'2  = 

7. 9967E-01 

WnT* 

6.91 13E-01 

DET*-2. 644OE-01 

P  1  * 

2. 2760E+O0 

P2* 

7.6 7396-01 

WnT  = 

7. 3398E-01 

DET  *~3.  1743E-01 

PI* 

2 . O643E+0O 

P2* 

7.  4656E-01 

WnT* 

3. 1SSIE-01 

D  E  T  * - 3 . 7269E-3 1 

PI* 

1 . 9108E  +  00 

P2* 

7. 3324E-01 

WnT* 

3. 7965E-0 1 

DET=-4. 2924E-01 

PI* 

1 . 7939E+00 

P2* 

7. 2617E-01 

WnT* 

9. 424SE-01 

DET*-4. 3615E-01 

PI* 

1 . 7  1  7  6  E  +  0  0 

P2* 

7.24016-O1 

WnT* 

1 . 0053E+00 

DET*-5. 42516-01 

PI* 

1 . 6390E+O8 

P2* 

7 . 2572E-0 1 

WnT  * 

1 . 0631 E+00 

D6T*-5 . 9746E-0 1 

PI* 

1.61 77E+0O 

P2* 

7. 305OE-01 

WnT* 

1 . 131OE+00 

DET*-6. 5820E-01 

PI* 

1 . 5S95E+80 

P2  = 

7. 3767E-0 1 

WnT* 

1 . 1 933E+O0 

D£T*-7. 0000E-0  1 

PI* 

1 . 3713E  +  00 

P2* 

7. 4672E-01 

WnT* 

1 . 25666+00 

DET*-7. 4622E-01 

PI* 

1.561 OE+00 

P2* 

7 . 5720E-0 1 

WnT  = 

1 . 3 1 95E+O0 

DET*-7.3830E-01 

PI* 

1 . 5566E+80 

P2* 

7. 6S75E-01 

WnT* 

1 . 3823E+00 

DET=-S.2S77E-01 

PI* 

1 . 5569E+00 

P2* 

7.31 06E-O 1 

WnT* 

1 . 443 l £+00 

DET*-$. 5827E-01 

PI* 

1 . 5607E+0© 

P2* 

7. 9336E-01 

WnT* 

1 . 5O30E+00 

DET—3.3349E-01 

PI* 

1 . 36726+00 

P2* 

3.06956-01 

WnT  * 

1 . 570SE+00 

DET*-9. 07276-01 

PI* 

1 . 5757E+00 

P2* 

S.2013E-01 

WnT* 

1 . 6336E+00 

DET--9. 2343E-01 

PI* 

1 . 5357E  +  00 

P2* 

8. 3324E-01 

WnT* 

1 . 69656+00 

DET*-9.3410E-01 

PI* 

1 . 3966E+00 

P2* 

3.461 6E-0 1 

WnT* 

1 . 73936+00 

DET*-9 . 39 1 9E-0 1 

PI* 

1 . 6OS2E+00 

P2* 

3.53766-01 

WnT* 

1 . 322  IE +00 

DET  *-9 . 33366-01 

PI* 

1 . 620 1 E  +  00 

P2* 

8. 7096E-01 

WnT* 

1 . 3S50E+00 

DET  =  -9 . 3323E-0 1 

PI* 

1 . 6321 E+00 

P2  = 

3 . 326SE-0 1 

WnT  * 

1 . 9  4  7  3  E + 0  0 

DET* -9 . 2269E-0  1 

PI* 

1 . 6440E+0O 

P2* 

3.93346-01 

WnT* 

2. O106E+00 

DET  *-9. 0733E-01 

PI* 

1 . 6556E+00 

P2* 

9. 0440E-01 

WnT* 

2. 0735E+00 

D  E  T  * - 3 . 3765E-01 

PI* 

1 . 6669E+08 

P2* 

9.  2  432E-01 

WnT* 

2. 1 363E+O0 

DET*-8. 6383E-01 

PI* 

1.6777E+O0 

P2* 

9.23376-01 

WnT* 

2.  1991E  +  00 

DET  *-3. 36376-01 

PI* 

1 . 6380E  +  O0 

P2* 

9. 32116-01 

WnT* 

2. 26196+00 

DET=-3 . O333E-0 l 

PI* 

l . 6976E  +  O0 

P2* 

9. 3995E-01 

WnT  = 

2. 32486+00 

DET*-7. 71906-01 

Pi* 

1 . 7066E+0O 

P  2* 

9. 470SE-0 1 

WnT* 

2 . 3376E+00 

DET--7. 3572E-01 

PI* 

1 . 7 1 43E+O0 

P2* 

9. 5348E-0 1 

WnT* 

2. 43O4E+00 

DET*-6. 9746E-01 

PI* 

1 . 72246+00 

P2* 

9.3913E-01 

WnT* 

2.31 33E+00 

DET  *-6 . 3750E-0  1 

PI* 

t .  7293E+00 

P2» 

9. 6417E-01 

WnT* 

2. 3761E+Q0 

DET*-6. 1623E-01 

PI* 

1 . 7354E+00 

P2* 

9 . 634SE-0 1 

WnT* 

2.6339E+00 

DET*-5. 7407E-01 

PI* 

1 . 74096+00 

P2* 

9.721 2E-0 1 

WnT* 

2. 7013E+00 

DET*-5. 3133E-01 

PI* 

1 . 7436E+O0 

P2* 

9. 751  IE-01 

WnT* 

2. 7646E+O0 

DET*-4. 3335E-01 

PI* 

1 . 7497E  +  00 

P2* 

9. 774SE-01 

WnT* 

2. 3274E+0O 

DET=-4. 4347E-0 1 

PI* 

1 . 7331  E+00 

P2* 

9. 7926E-01 

WnT* 

2. 3903E+00 

DET  *-4 . 0296E-0 1 

Pi* 

1 . 7359E+00 

P2* 

9.30476-01 

WnT* 

2 . 933 16+00 

DET*-3.61 10E-01 

PI* 

1 . 7531 E  +  00 

P2* 

9.31 16E-01 

WnT* 

3.01 59E+00 

DET*-3. 2013E-O1 

PI* 

1 . 7397E+0O 

Pit!  * 

9.31346-01 

WnT* 

3 . O733E+00 

DET  *-2 . 38266-0 1 

PI* 

1 . 7  6O3E  +  00 

P2  = 

9.31066-01 

WnT  * 

3. 1416E+00 

DET*- L . 4167E-01 

PI* 

1 . 7614E  +  00 

P2* 

9.30356-01 

2E  TR 

*  7.O700E-01 

INITIAL  ONEGA  «  0. 

0000E+00 

FINAL  OMEGA  *  3.1416E+00 

NUMBER  OF  STEPS* 

50 

W  n  T  ■ 

:  1.0000E-23 

THE  DETERMINANT  OF 

THE  COEFFICIENT 

MATRIX 

IS  EQUAL  TO  ZERO 

MATRIX  AB  » 

O. 0000E+00 

0. 0000E+00 

0 . 0000E+00 

1.41 44E-25 

MATRIX  Fp  * 

-3. 0000E-01 

-5. 6577E-26 

WnT* 

6.2332E-02 

DET  *-3 . 2092E-04 

PI* 

2 . 0273E+O2 

P2* 

9 . O067E+00 

WnT  * 

1 . 2566E-01 

BET  —  2. 3468E-03 

PI* 

5 . 0  3  3 1 E  +  0 1 

P2* 

4.3141 E+00 

WnT* 

1 . S350E-0 1 

DET  *-7. 235 1 E-03 

PI* 

2. 2766E  +  01 

P2  = 

3. 0250E+00 

WnT* 

2.51 33E-0 1 

DET  =  - 1 . 5656E-02 

PI* 

1 . 2991E  +  01 

P  2* 

2. 23S3E+00 

WnT* 

3.  1416E-01 

DET*-2. 7895E-02 

PI* 

3. 504SE  +  00 

P2* 

l .  3547E  +  00 

WnT* 

3. 7699E-01 

DET=-4. 3945E-02 

PI* 

6. O977E+00 

P2* 

1 . 5725E+00 

WnT* 

4 . 39S2E-0 1 

DET=-6 . 3573E-02 

PI* 

4. 6710E+00 

P2* 

1 . 3776E+O0 

WnT* 

5. 02S5E-01 

DET--3.6406E-02 

PI* 

3. 7S55E+00 

P2* 

1 . 2376E+00 

WnT* 

5. 6549E-01 

DET  =  - 1 . 1 194E-01 

PI* 

3.  1S19E+00 

P2* 

1 . 1 342E+00 

WnT* 

6 . 2832E-0 1 

DET*- 1 . 3961 E-0 1 

PI* 

2. 7443E  +  &yi 

P2* 

1 . 0566E+00 

WnT* 

6. 91 15E-01 

DET=-i . 6835E-0 1 

PI* 

2 . 4  4  3  6  E + O  0 

P2* 

9 . 9772E-0 1 

WnT* 

7. 3398E-01 

DET  =  - 1 . 9905E-O 1 

Pi* 

2. 234  1E  +  0O 

P2  = 

9 . 5233E-0 1 

WnT  = 

3. 1 631 E-0 1 

DET=-2. 2964E-01 

PI* 

2 . O765E+00 

P2  = 

9.  1S67E-01 

WnT  * 

3. 7965E-01 

DET  *-2 . 300SE-0 1 

PI* 

1 . 9594E  +  00 

P2* 

3. 9273E-01 

WnT* 

9. 4243E-01 

DET  *-2 . 8987E-0 1 

PI* 

1 . 8719E+00 

P2* 

3. 7341E-01 

WnT* 

1 . O053E+00 

DET  *-3 . 1357E-01 

PI* 

1 . 3085E  +  00 

P2* 

3. 3920E-0 1 

WnT* 

1 . 0631 E+00 

DET*-3. 4577E-01 

PI- 

1 . 7575E+O0 

P2* 

8.491  IE-01 

WnT* 

1 . 1310E+00 

DET*-3. 71 1SE-01 

PI* 

1.721 1 E+00 

P2* 

8 . 423 1 E-0 1 

WnT* 

1 . 1 933E+00 

DET=-3. 9445E-01 

PI* 

l .  S943E+00 

P2* 

8.3816E-01 

WnT* 

1 . 2366E+00 

DET  *-4 . 1543E-01 

PI* 

1 .3750E+00 

P2* 

8.361  IE-01 

WnT* 

1 . 3195E+00 

DET  *-4 , 3394E-0 1 

PI* 

1 . 6S14E  +  00 

P2* 

8 . 3576E-0 1 

WnT* 

1 . 3323E+00 

DET=-4. 4937E-01 

PI* 

1 . 6523E+00 

P2* 

8. 3674E-81 

WnT  = 

1 . 4451E  +  00 

DET*-4. 6315E-01 

PI* 

1 . 3 437E+O0 

P2  = 

3. 3S76E-01 

WnT* 

1 . 5080E+00 

DET=-4. 7377E-01 

PI* 

1 . 6437E  +  00 

P2* 

3.41 5SE-0 1 

WnT* 

1 . 37Q3E+00 

DET=-4. 3174E-01 

PI* 

1 . S429E+00 

P2* 

3. 4501E-01 

WnT  = 

1 . 6336E+00 

DET*-4. 3712E-01 

PI* 

1 . S43SE+00 

P2* 

3 . 483SE-8 1 

WnT* 

1 . 6965E+00 

DET  *-4 . S999E-0 1 

PI* 

1 . 3455E  +  00 

P2* 

8 . 5306E-0 1 

WnT  = 

1 . 7593E+00 

DET  =  -4 . 9045E-O 1 

PI* 

1 . S4S3E+00 

P2* 

8. 5742E-0 1 

WnT* 

1 . 322  IE +  00 

DET  =  -4 . 8364E-0 1 

P  l  * 

1 . 65 1 3E  +  0O 

P2* 

8.61 83E-0 1 

WnT* 

1 . 3350E+00 

DET*-4. 8469E-01 

PI* 

1 . 6557E+00 

P2* 

8 . 6636E-0 1 

WnT* 

1 . 9478E+80 

DET=-4. 7375E-01 

PI* 

1 . 6593E+0O 

P2  = 

3. 7079E-01 

WnT* 

2.01 06E+00 

D  E  T  =  - 4 . 7099E-01 

PI* 

1 . 664 1 E  +  00 

P2* 

8. 7512E-01 

WnT  = 

2 . 0735E+00 

DET--4.6159E-01 

PI* 

1 . 6634E  +  00 

P2  = 

3. 7932E-0 1 

WnT  = 

2. 1 363E+00 

DET=-4. 507OE-0 1 

PI* 

1 . 6727E+0O 

P2* 

3. 8334E-01 

WnT* 

2.  1991E+00 

DET*-4. 3350E-01 

PI* 

1 . 6769E+00 

P2* 

8. 8717E-01 

WnT* 

2. 2S19E+00 

DET*-4.2516E-01 

PI* 

1 . 6310E+00 

P2* 

3 . 907SE-0 1 

WnT* 

2 . 3248E+00 

DET  =  -4.  1034E-01 

PI* 

1 . 6S49E+00 

P2* 

3. 9416E-01 

WnT  * 

2. 337SE+00 

DET  »-3 . 957 1 E-0 1 

PI* 

1 . 68S5E  +  00 

P2* 

8. 9731E-01 

WnT  * 

2. 4504E+00 

DET  =  -3 . 799 1 E-0  1 

PI* 

1 . 6919E  +  00 

P2* 

9 . 0022E-O 1 

WnT  * 

2.51 33E+00 

DET  *-3. 6359E-0  1 

PI* 

1.695 1 E+00 

P2* 

9 . 0288E-0 1 

WnT  * 

2. 57S1E+00 

DET--3. 4689E-01 

PI* 

1 . 6980E+00 

P2  = 

9 . 0530E-O 1 

WnT  * 

2. 6389E+00 

DET  *-3 . 2994E-0 1 

PI* 

1 . 70O7E+00 

P2* 

9.  O743E-01 

WnT* 

2.701 3E+0O 

DET  *-3. 1 235E-0 1 

PI* 

1 . 7032E  +  O0 

P2* 

9. 0944E-O1 

WnT* 

2 . 7646E+00 

DET—2. 9574E-01 

PI* 

t . 7034E+00 

P2* 

9.  1 1 17E-01 

Wn  T  * 

2. 3274E+00 

DET*-2. 7339E-01 

PI* 

1 . 7073E+0O 

P2* 

9.  1 26SE-0 1 

WnT* 

2. 3903E+00 

DET  *-2 . 6 1 3 1 E-0 1 

PI* 

1 . 7090E  +  00 

P2* 

9.  1399E-01 

WnT *  2.9531E+00  DET*-2. 4516E-01  Pi  *  1.7105 E+09  P2=  9.1510E-O1 
WnT*  3.91 59E+00  DET--2. 2333E-01  *  Pi*  1.7113E+00  P2=  3.  1603E-O1 
W n T *  3.073SE+0O  DET*-2. 1287E-01  Pi*  1.7129E+09  P2*  9.1679E-01 
WnT*  3.141  6E+O0  DET *- 1 .  9733E-0 1  PI*  1 . 7 1  39E+00  P2=  9.1739E-01 


ZETA  *  9.9999E-01S1 

INITIAL  OMEGA  *  0.0000E+00  FINAL  OMEGA  =  3.1416E+00 

NUMBER  OF  STEPS*  50 
WnT  *  1 . 0000E-25 

THE  DETERMINANT  OF  THE  COEFFICIENT  MATRIX  IS  EOUAL  TO  ZERO 
MATRIX  AB  * 

0.0000E+00  0.0000E+00 

0.0000E+00  S.9442E-2S 

MATRIX  Fp  * 

-3. 0000E-01  -3 . 5777E-28 


WnT  = 

S.2S32E-02 

DET  *- 1 . 9573E-06 

PI* 

2 .  O  9  O  £  +  O 

P2* 

1 . 273SE+0 1 

WnT* 

l .  2566E-0 1 

DET  *- 1 . 3823E-05 

PI* 

5.  1 023E+0 1 

P2* 

6. 331 5E +O0 

WnT* 

1.3350E-O1 

DET*-4. 121  IE-05 

Pi* 

2. 2973E+01 

P2* 

4. 2714E+00 

WnT* 

2.51 33E-0 1 

DET  *-3. 6348E-05 

PI* 

1 . 3213E+01 

P2* 

3. 2245E+O0 

WnT* 

3.  1416E-01 

DET*- 1 . 49 1 7E-04 

PI* 

3. 7422E+00 

P2* 

2. 6037E+00 

WnT* 

3. 7699E-01 

DE  T  * - 2 . 28 1 5E - 8  4 

PI* 

6 . 341 6E+00 

P2  = 

2.  1 964E  +  O0 

WnT* 

4 . 3  9  8  2  E  -  0 1 

DET*-3 . 2088E-04 

PI* 

4.91 79E+O0 

P2* 

1 . 91 14E+0O 

WnT* 

5.0265E-O1 

DET*-4.2450E-04 

PI* 

4.01 26E+00 

P2  = 

1 . 7O29E+00 

WnT* 

5. €5496-01 

DET*-5. 36O2E-04 

F'l* 

3. 4067E+00 

P2* 

1 . 5455E+00 

WnT  = 

6.2S32E-01 

DET* -6. 5250E-04 

Pi* 

2.9353E+O0 

P2* 

1 . 4237E  +  O0 

WnT  * 

6.911 5E-0 1 

DET*-7.7U9E-04 

PI* 

2 . 683 1 E+00 

P2* 

1 . 3277E+0O 

WnT* 

7. 5393E-01 

0ET*-3. 3963E-04 

PI* 

2.4612E+00 

P2* 

1.251 0E+00 

WnT  = 

3. 1631E-01 

DET*- l . 0057E-03 

PI* 

2 , 295 1 E+00 

P2* 

1 . 1 390E+00 

WnT  * 

3. 7965E-0 1 

DET  *- 1 . 1175E-03 

PI* 

2. 1 687E+00 

P2* 

1 .  1  385E  +  00 

WnT* 

9.4243E-01 

DET*-1 .2236E-03 

PI* 

2.071 2E+00 

P2* 

1 . 0970E+00 

WnT* 

1 . 0053E+O0 

DET  *- 1 . 3229E-03 

PI* 

1 . 9953E  +  0O 

P2* 

l . 0627E+00 

WnT* 

1 . 063 1 E  +  00 

DET=-1.4144E-03 

PI* 

1 . 9355E+0O 

P2* 

l . 0342E+O0 

WnT* 

1 .  1310E+00 

DET*- 1 . 4975E-03 

PI* 

1 . 3330E+00 

P2  = 

l . 0105E+O0 

WnT* 

1 . 1 933E+00 

DET*- 1 . 5713E-03 

PI* 

1.3502E+00 

P2  = 

9. 9074E-O1 

li 

1 . 2566E+0O 

DET*- 1 . 6372E-03 

PI* 

1 . 3193E+0O 

P2  = 

9. 7416E-01 

WnT* 

1 . 3 1 95E+0O 

DET*-1 .6935E-03 

PI* 

1 . 7953E+00 

P  2* 

9. 6025E-01 

WnT* 

1 . 3323E+0O 

DET*- 1 . 7410E-03 

PI* 

1 . 7754E+0O 

P2  = 

9. 4859E-01 

WnT  = 

1 . 4451 E+00 

DET  *- 1 . 7797E-03 

PI* 

1 . 7593E+0O 

P2* 

9. 3S79E-01 

WnT* 

l . 503OE+00 

DET*-1 . 8101E-03 

PI* 

1 . 7462E  +  00 

P2  = 

9. 3057E-0 1 

WnT* 

l . 5703E+O0 

DET=- 1 . S326E-03 

PI* 

1 . 7356E  +  0O 

P2  = 

9. 2367E-01 

WnT* 

1 . 6336E+0O 

DET*- 1 . 3475E-03 

PI* 

1 . 726 9 E+00 

P2  = 

9.  1  789E-0 1 

WnT  = 

1 . 6965E  +  0O 

DET*-1 .8553E-03 

PI* 

1 . 7198E+00 

P2  = 

9.  1  306E-0 1 

WnT* 

1 . 7593E+O0 

DET*- 1 . 3566 E -03 

PI* 

1 . 714OE+00 

P2* 

9. O903E-01 

WnT* 

1 . 822 1 E+00 

DET*- 1.351 9E-03 

Pi* 

l  ..7093E+00 

P2* 

9. 056SE-01 

WnT* 

1 . 335OE+00 

DET  *-1.341 6E-03 

PI* 

1 . 7055E+0O 

P2* 

9. O290E-01 

WnT* 

1.9473E+O0 

DET*-1.3262E-03 

PI* 

1. 7024E+0O 

P2* 

9. 0062E-O 1 

WnT  * 

2. 01O6E+00 

DET*-l . 3064E-03 

Pl- 

1 . 6999E+0O 

P2* 

3. 9375E-01 

WnT* 

2. 0735E+0O 

DET  *- 1 . 7S25E-03 

Pl* 

1 . 69SOE+00 

P2* 

8. 9723E-01 

WnT* 

2 . 1 363E+00 

DET  *- 1 , 7553E-03 

Pi* 

1 . 6964E+00 

P2* 

8.9600E-O1 

WnT* 

2. 1991 E+00 

DET*- 1 . 7244E-03 

pi* 

1.6952E+00 

P2* 

3. 9503E-01 

WnT  * 

2. 2619E+00 

DET  *- 1 . 6 9 1 0  E - O  3 

pi* 

1 . 6942E+0O 

P2* 

3. 9423E-01 

WnT* 

2. 3243E+00 

DET*-t . 6554E-03 

pi* 

1 . 6935E+0C 

P2* 

8. 9370E-0 1 

WnT* 

2 . 3876E  +0O 

DET  =  - 1 . 6 1 77E-03 

Pi* 

1 . 6930E+0O 

P2* 

3. 9323E-01 

WnT  * 

2. 4504E+30 

DET*- 1 . 5784E-03 

pi* 

t . 6927E+00 

P2- 

3 . 9293E-0 1 

WnT  = 

2.51 33E+00 

DET  =  -  1 . 

5379E-03 

pi  = 

1 . ©924E+00 

P2* 

8.9279E-01 

WnT  = 

2.57S1E+00 

DET=- 1 . 

4962E-03 

pi* 

1 . G923E+00 

P2* 

3. 92E9E-01 

WnT  = 

2 . S389E+00 

DET  — 1 . 

4  5  3  9  E  -  0  3 

pi* 

1 . 6923E+00 

P2= 

8. 9266E-01 

WnT  = 

2.701 3E+00 

DET  — 1 . 

41  10E-03 

pi* 

1 . G923E+00 

P2* 

3. 9269E-01 

WnT* 

2 . 7S46E+00 

DET— 1 . 

3673E-03 

pi* 

1 . G924E+00 

P2* 

8. 9277E-01 

WnT  = 

2. 3274E+00 

DET  — 1. 

3244E-03 

pi* 

1 . 692GE+00 

P2  = 

3. 9239E-01 

WnT  = 

2.S903E+00 

DET— 1. 

2S12E-03 

pi* 

1 . 6927E+00 

P2* 

8 . 9305E-0 1 

WnT  = 

2. 9531E+00 

DET—  1 . 

23S1E-03 

pi* 

1 . G929E+O0 

P2* 

3 . 9322E-0 1 

WnT  = 

3. 0159E+00 

DET  —  1  . 

1954E-03 

pi* 

1 . 6931E+00 

P2* 

3. 9342E-0 1 

WnT  * 

3 . 0733E+00 

DET  —  1 . 

1532E-03 

pi* 

1 . S934E  +  00 

P2* 

8 . 9363E-0 1 

WnT  = 

3.  1416E+00 

DET  —  1 . 

1 1  loE-03 

pi* 

1 . G936E+00 

P2* 

3 . 9335E-0 1 

Now  a  tentative  design  point 

may 

be  selected. 

Let  us  suppose 

that 

,  for  the  first 

desi gn 

iteration. 

the 

control  system  designer  selects 

a  location  of  5  =  0 

.707  and 

w  T  3  ir/2 

for 

the  placement  of  the  poles 

.  This 

corresponds  to  a  point  in  the  parameter  space 

at  ap  *  0.845 

,  a.  3  0. 

2,  and 

aD  = 

1  .643.  Notice 

that  this  establishes 

the 

value  for  the 

real  root,  i.e.. 

6  = 

0.02.  If  other 

values  of  5  were 

desi red, 

the  value  of 

C  3  0.707 

could 

be  maintained  and  the  value 

of  id  T  varied 

n 

or 

vice  versa. 

These  values  for 

root  locations  may  be  independently  checked  by  using  the  programmed  poly¬ 
nomial  root  finder  termed,  by  Hewlett-Packard,  "Siljak's  Method"  (Appendix 
F). 


COEFFICIENTS  OF  THE  CHARACTER  I  ST  I C  EQUATION  WITH 


pi*  i.< 

5430E+00 

P2*  8. 4500E-01 

c  <  0  >  + 

C  < 1 >  *S  + 

C<2>*$/'2 

+  C  <  3  >*3''  3  +  ... 

c  <.  0  > » ■ 

■2.00O000E 

-03 

O 

*-• 

V' 

II 

1 . 143033E 

-01 

c«:  2)*- 

-3. 120000E 

-01 

c  <  3  )  * 

1 . 000000E 

+  00 

ROOTS: 

REAL  IMAGINARY  DAMPING  RATIO 


1 . 467913E-01 
I . 34 1 74OE-02 
1 . 4S7913E-01 


-2. 950344E-0 1 
1.608759E-13 
2.930344E-01 


7.0741 1  IE-01 
1 , 000O00E+00 
7.0741  1  IE-01 


WnT 


1 . 569207E+00 
3 . 994480E+00 
1 . 5S9207E+00 
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Since  the  c-contours  by  interpolation  span  the  entire  stable 
region,  there  must  exist  (in  that  region)  only  one  real  root  at  each  point 
If  real  root  locations  in  the  z-domain  are  designated  by  the  symbol  6,  the 
6-contours  (straight  lines)  are  solved  for  values  of  6  =  0.8,  0.7,  0.5, 
0.2,  -0.1,  -0.5,  and  -0.9  (notice  of  course,  that  the  z  -  -1  stability 
boundary  is  synonymous  with  the  5  *  -1  contour). 

FOR  REAL  ROOTS  OF  THE  CHARACTERISTIC  EQUATION 

DELTA  =  8.0000E-O1 

P2  =  1.111  IE-01  *P1  +  1 . 7778E+0O 

THE  AXES  INTERCEPTS  ARE  AS  FOLLOWS: 

AT  P1=0  ,  P2*  1.7778E+00 

AND  AT  P2-0  ,  P 1 *- 1 . SO00E+O1 


FOR  REAL  ROOTS  OF  THE  CHARACTERISTIC  EQUATION 

DELTA  =  7.000OE-01 

P2  =  1 . 7647E-0 1  *P1  +  1 . 0804E+0O 

THE  AXES  INTERCEPTS  ARE  AS  FOLLOWS: 

AT  Pl*e  ,  P2*  1 . 0304E+0O 
AND  AT  P 2*0  ,  P 1 *-€ . 1 222E+00 


FOR  REAL  ROOTS  OF  THE  CHARACTERISTIC  EQUATION 

DELTA  *  5.OO00E-01 

P2*  3. 3333E-01  *P1  +  4.3333E-01 

THE  AXES  INTERCEPTS  ARE  AS  FOLLOWS: 

AT  P1»0  ,  P2=  4 . 3333E-0 1 
AND  AT  P2-0  ,  P 1 *- 1 . 3O00E+00 


FOR  REAL  ROOTS  OF  THE  CHARACTERISTIC  EQUATION 

DELTA  »  2.0000E-01 

P2-  <5. 66S7E-01  *P1  ♦  -2.3333E-01 

THE  AXES  INTERCEPTS  ARE  AS  FOLLOWS: 

AT  P 1 *0  ,  P2*-2 . 3333E-0 1 
AND  AT  P2*0  ,  PI*  3.500OE-01 
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FOR  REAL  ROOTS  OF  THE  CHARACTERISTIC  EQUATION 

DELTA  *-1.0000E-0t 

P2»  1.2222E+00  *Pl  *  -1.1S0SE+00 

THE  AXES  INTERCEPTS  ARE  AS  FOLLOWS: 

AT  P1*0  ,  P2  =  - 1 . 1 80SE+00 


AND  AT  P2-0  ,  PI*  9.SS12E-01 


FOR  REAL  ROOTS  OF  THE  CHARACTERISTIC  EQUATION 

DELTA  *-5 . 0000E-0 1 

P2*  3.O000E+00  *Pl  +  -4. 4333E+00 

THE  AXES  INTERCEPTS  ARE  AS  FOLLOWS: 

AT  Pl30  ,  P2=-4. 4333E+00 
AND  AT  P2*0  ,  PI*  1 . 4778E+00 


FOR  REAL  ROOTS  OF  THE  CHARACTERISTIC  EQUATION 

DELTA  *-9. 0000E-01 

P2*  l . 9000E+0 1  *Pl  +  -3 . S089E+0 1 

THE  AXES  INTERCEPTS  ARE  AS  FOLLOWS: 

AT  P 1*0  ,  P2--3. 6089E+01 
AND  AT  P2*0  ,  PI*  I . 8994E+O0 


4.  Application  of  Cross-Multiplication  Technique 

Suppose  now  the  control  system  designer  wishes  to  assess  the 
system  response,  using  the  numerical  values  for  ap,  a^,  and  a^.  Since  he 
already  has  the  closed-loop  transfer  function  of  Equation  (20),  he  may  go 
directly  to  the  Cross-Multiplication  method.  The  digital  program  for  this 
method  is  shown  in  Appendix  G.  He  loads  in  the  order  of  numerator  (2)  and 
denominator  (3)  and  selects  a  value  for  sampling  interval  (here,  T  *  Is). 
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ORDER  OF  NUMERATOR  IS  2  ORDER  OF  DENOMINATOR  IS  3 

SAMPLE  INTERVAL  =  1 


CALCULATIONS  MADE  ACCORDING  TO  FOLLOWING  EQUATION  DEFINITION 
C  C  Z  >  '•'  R  <  2  >  =  NUMERATOR-DENOMINATOR 

NUMERATOR  =  AO3>*Z-'0  ♦  fl<l>»Z-l  +  A<.2>*ZA2  .  .  .  A < N >  * Z N 
DENOMINATOR  =  B<:O>*Z'-0  ♦  EKl)*Z'l  +  B<2>#2*2  .  •  .  £•  D>*Z-  D 


Next,  the  numerical  values  for  the  coefficients  in  the  numerator  and 
denominator  are  loaded. 

NUMERATOR  COEFFICIENTS  H ■ a > , A <  1 > . A t  2 • . . . . A'  N> 

. 998  -2.886  2.633 


DENOMINATOR  COEFFICIENTS  B<0> ,B< 1> , B<2)  ,  . . .  B<N> 
-.002  .114  -.312  1 


Then,  the  number  of  sampling  instants  to  be  used  is  selected.  In  this 
case,  time  was  varied  from  zero  to  10T.  The  test  input  selected  for  ec 
(here  designated  "INPUT")  is  a  unit  step  function.  The  output  e  is  printed 
out  as  "RESPONSE."  It  is  plotted  as  Figure  3.  The  R’s  indicate  values 
for  the  input  at  sampl ing  instants,  and  the  asterisks  indicate  value  for 
the  response  at  the  same  instants. 

NUMBER  OF  STEPS  FOR  RESPONSE  =  10 
THE  INPUT  IS  A  UNIT  STEP  FUNCTION 


ec 

0 

K 

TINE 

INPUT 

RESPONSE 

0 

0. 000E+00 

1 . 000E  +  00 

0 . 000E+00 

1 

1 . 000E+00 

1 . 000E+00 

2. 688E+0O 

2 

2 . 000E+0O 

1 . 000E+00 

6. 407E-01 

3 

3.  O00E  +  00 

1 . O00E+00 

6.935E-01 

4 

4 . 000E+00 

1 . 000E+00 

9. 487E-01 

5 

5.000E+00 

1 . 000E  +  00 

1 . 018E+00 

6 

6. 000E+00 

1 . 000E+00 

1.01 1E+00 

7 

7. 000E+00 

1.000E+O0 

1.001E+00 

3 

3 . 000E+00 

1 . 000E+00 

9 . 992E-0 1 

9 

9. 000E+00 

1 . 000E+00 

9. 996E-01 

10 

1 . 000E+01 

1 . 000E+00 

1 . 000E+80 

At  this  point  the  designer  may  wish  to  iterate  upon  his  selected  design 
point,  with  his  frequency  and  time-domain  information  now  in  hand. 
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SECTION  IV.  DIGITAL  SYSTEM  SEMINAR  (Task  V) 


t 


! 


>  i  I 


n 


During  the  December  1979  -  September  1980  portion  of  the  contracted 
period,  a  Digital  Design  Seminar  was  conducted  for  selected  members  of  the 
Guidance  and  Control  Directorate.  The  text.  Digital  Control  Systems,  by 
8.  C.  Kuo  was  provided  the  participants  by  the  G&C  Directorate. A  detailed 
outline  of  the  topics  covered  follows. 

A.  Introduction 

1.  Objective  of  seminar 

2.  Philosophy  of  seminar 

3.  Why  digital  control? 

4.  Problem  definition 

5.  Description  of  the  sampling  operation 

6.  The  Z- transform 

7.  Examples  of  analysis  of  digital  systems  by  "cook  book  recipes" 

a.  SAM 

b.  Stability 

c.  Response  by  Cross-Multiplication 

B.  The  sampling  process 

1.  Motivation 

2.  Mathematical  description  by  Pulse  Amplitude  Modulation 

a.  Time-domain  approach 

b.  Frequency  spectrum  approach 

c.  Complex  convolution  approaches 

C.  The  Sampling  Theorem 

1.  Limits  on  sampling  frequency 

2.  Shannon's  Theorem 

3.  Fogel 

D.  The  "Ideal  Sampler":  Impulse  sampling  approximation 

1.  Motivation 

2.  Mathematical  descriptions 

a.  Time-domain 

b.  Frequency  domain  using  Fourier  analysis 


c.  Frequency  domain  using  Laplace  analysis 

(1)  Series  forms 

(2)  Closed  form 

3.  Zero  Order  Hold 

4.  S-domain  properties 

E.  Z-transform  analysis  of  linear  sampled-data  systems 

1.  Introduction 

a.  Series  form 

b.  Closed  form 

c.  Advantages  versus  disadvantages 

2.  Evaluation  of  z-transforms 

3.  Mapping  the  s-plane  onto  the  z-plane 

4.  Inverse  z-transform 

a.  Inversion 

b.  Methods 

(1)  Partial  Fraction 

(2)  Power  Series 

(3)  Inversion  Formula 

(4)  Cross-Multiplication 

5.  Theorems  of  the  z-transform 

a.  Linearity 

b.  Real  translation 

c.  Complex  translation 

d.  Final  value 

e.  Initial  value 

f.  Partial  differentiation 

6.  Pulse  transfer  functions 

7.  Limitations  of  z-transform  method 

8.  Response  of  open-loop  sampled-data  systems  between  sampling 

instants 

a.  Discussion 

b.  Submultiple  sampling  method 

c.  Modified  z-transform 
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d.  The  inverse  modified  z-transform 

(1)  Introduction 

(2)  Power  series  expansion  series  method 

(3)  Inversion  formula 

(4)  Cross-multiplication  method 
9.  Modified  z-transform  theorems 

a.  Real  translation 

b.  Complex  translation 

c.  Initial  value  theorem 

d.  Final  value  theorem 

F.  Data  reconstruction 

1.  Introduction 

2.  Data  reconstruction  by  polynomial  extrapolation 

3.  Zero-order  hold 

4.  First-order  hold 

5.  Other  holds 

6.  Summary 

6.  Sampled-data  systems 

1.  Introduction 

2.  Block  diagram  analysis  and  transfer  functions  of  closed-loop 
sampled-data  systems 

a.  Introduction 

b.  Transfer  functions  and  output  response  of  open-loop 
systems  with  cascaded  elements 

c.  Closed-loop  sampled-data  systems 

d.  Summary  of  some  appropriate  observations 

e.  Additional  comments:  parallel  paths 

f.  Motivation  for  Signal  Flow  Graph  (SFG)/SAM  methods 

3.  Signal  Flow  Graphs  (SFG) 

a.  Introduction 

b.  The  sampled  SFG  method 

c.  "Formal"  SFG  method 

d.  Modified  SFG  method 

e.  SAM:  "Systematic  Analytical  Method" 
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4.  Modified  z-trartsforms  of  outputs  of  closed-loop  sampled- 
data  systems 

H.  System  response 

1.  Introduction 

2.  Comparison  of  transient  response  and  sampled  systems 

3.  Stability  of  sampled-data  systems 

a.  Introduction 

b.  Hidden  oscillations 

c.  Review  of  continuous-data  system  stability 

d.  Sampled-data  system  stability 

e.  Characteristic  equation 

f.  Schur-Cohn  stability  criterion 

g.  Jury's  stability  test 

(1)  Description 

(2)  Strengths  and  weaknesses 

(3)  Singular  cases 

h.  Raible's  tabular  form 

(1)  Description 

(2)  Strengths  and  weaknesses  — 

(3)  Singular  cases 

i.  Extension  of  Routh-Hurwitz  Criterion  to  sampled-data 

systems 

(1)  In  Laplace  (s)  domain 

(2)  Bilinear  transformations 
•r 

(a)  w-transform 

(b)  r-transform 

(c)  Relation  between  real  frequency  and  warped  w-- 
or  r-domain  frequencies 

4.  Effect  of  pole-zero  configurations  (of  starred  and  z  functions) 
upon  system  transient  response 

5.  Effect  of  pole-zero  configurations  on  maximum  overshoot  and 
peak  time  of  transient  response  of  sampled-data  systems 

6.  Root  Locus  method  for  sampled-data  systems 


a.  In  Laplace  (s)  domain 

b.  In  z-domain 

c.  Rules  of  construction 

d.  Stability  determination 

e.  Transient  response  and  relative  stability 

f.  Discussion  of  method 

7.  System  response  between  sampling  instants 

a.  Introduction 

b.  Submultiple  sampling  method 

(1)  Formal  method 

(2)  Cross-multiplication  method 

c.  Modified  z-transform  method 

(1)  Power  series  method 

(2)  Inversion  Formular  method 

(3)  Cross-multiplication  method 

d.  Determination  of  hidden  instabilities  and  oscillations 

8.  Steady-state  response  and  steady-state  error  analysis 

I.  Frequency  response  of  sampled-data  systems 

1.  Introduction 

2.  Nyquist  criterion  and  plot 

a.  2-  or  3- term  approximation  method  (s-domain) 

b.  z-transform  method 

c.  Bilinear  transform  method(s) 

3.  Bode  diagram  (w  and  r  domains) 

4.  Gain  phase  plots  (z-domain) 

5.  Parameter  Space  method  (z-domain) 

a.  Stability  determination 

b.  Transient  response  determination 

c.  Computer  mechanization 

d.  Advantages  and  disadvantages 

J.  State  variable  (modern)  approaches 

1.  Introduction 

2.  Continuous-data  systems 
a.  State  equations 


b.  State  transition  equations 
3.  Digital  systems 

a.  Analogies  to  continuous-data  system  approaches 

b.  Discrete  state  equation 

c.  Discrete  state  transition  equations 

d.  Solution  of  time-invariant  discrete  state  equation  by 
z-transform  methods 

e.  Discussion  of  further  applications 


K.  Summary 
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SECTION  V.  PERSHING  II  DIGITAL  G&C  (Task  VI) 

The  applicability  of  applying  advanced  digital  technology  to  the 
PERSHING  II  G&C  system  was  investigated.  In  concert  with  Mr.  James  McLean 
(DRSMI-RGN) ,  the  digital  autopilot  was  selected  as  the  point  of  focus.  The 
results  of  these  investigations  were  provided  to  Mssrs.  Russell  Gambill  and 
McLean  in  the  form  of  oral  reports. 

Dr.  Seltzer  attended  and  participated  in  the  PERSHING  II  digital  auto¬ 
pilot  design  review  held  at  The  Martin  Company  in  Orlando  in  late  March 
1980.  Written  coimients  were  provided  to  Mr.  McLean  (Appendix  G).  An 
important  part  of  this  (Appendix  G)  document  is  the  concluding  portion  which 
is  comprised  of  six  recommendations,  repeated  below  for  continuity. 

1.  Provide  copies  of  the  proposed  presentation  material  to  all  reviewers 
approximately  one  week  before  the  next  technical  review. 

2.  Investigate  the  possible  existence  of  dynamically  significant  non- 
linearities.  Predict  their  dynamic  effect  by  mathematical  analysis  to  be 
verified  by  computer  simulation.  Where  this  is  deemed  impractical,  include 
in  the  hybrid  simulation  the  actual  hardware  that  gives  rise  to  the  non¬ 
linearity. 

3.  Investigate  analytically  the  dynamic  (and  accuracy)  effects  of 
digital  computer  quantization. 

4.  Investigate  the  sensitivity  of  the  control  system  to  variations  in 
numerical  values  of  bending  parameters.  View  the  5 %  accuracy  figure 
associated  with  the  NASTRAN  program  with  skepticism.  This  should  have  a 
strong  effect  on  the  planned  "notch  filters." 

5.  Investigate  the  sensitivity  of  the  control  system  to  realistic 
variations  in  numerical  values  of  system  parameters  involved  in  any  pole- 
zero  cancellation  techniques. 

6.  Initiate  action  to  decrease  computer  capacity  requirements  and/or 
strongly  limit  future  additional  demands  upon  the  onboard  computer. 
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SECTION  VI.  DIGITAL  FILTERS  (Task  VII) 


In  June  1980,  effort  was  begun  to  evaluate  the  consequences  of  micro¬ 
processor  implementation  of  digital  filters  and  the  attendant  control  logic. 
Three  key  issues  were  identified: 

1.  Accuracy  due  to  word  length, 

2.  Computational  capability,  and 

3.  Memory  requirements. 

It  was  determined  that  requirements  associated  with  the  identified  key 
issues  must  be  determined.  The  requirements  would  be  of  a  form  that  would 
then  enable  one  to  determine  hardware  and  control  implementation  require¬ 
ments.  Emphasis  will  be  from  a  digital  perspective,  as  opposed  to  attempt¬ 
ing  to  emulate  analog  systems  (as  is  apparently  the  current  thrust  in  con¬ 
temporary  1 i terature ) . 

An  investigation  of  implementation  techniques  and  required  analytical 
tools  required  for  digital  filtering  was  begun,  as  was  a  review  of  current 
microprocessor  state-of-the-art  and  hardware  characteristics.  This  was 
interrupted  by  the  arising  emphasis  on  ATBM  G&C  (Section  IIC). 
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SECTION  VII.  SUMMARY 

During  the  contracted  period  of  31  July  1979  through  30  September  1980, 
the  study  objective  set  forth  in  the  Contract  was  met,  i.e.,  to  perform 
research  necessary  to  develop  new  discrete  design  techniques,  discrete 
adaptations  of  control  theory,  and  in  particular  extended  Parameter  Space 
analysis  techniques  that  will  help  to  ensure  development  of  effective  G&C 
systems  for  future  weapon  systems.  The  study  effort  was  composed  of  seven 
tasks.  The  requirements  established  in  Task  I  were  performed  and  documented 
on  30  September  1979  in  a  separate  technical  report  entitled,  "A  Mathematical 
Means  of  Spreading  Clustered  System  Roots."  The  requirements  of  the  remain¬ 
ing  tasks  were  performed  and  are  documented  in  this  Final  Report. 

One  of  the  most  significant  technical  achievements  to  be  accomplished 
in  this  contracting  period  was  the  development  and  documentation  of  a  simple 
digital  system  design  tool.  This  tool  is  a  combination  of  the  z-domain 
Parameter  Space  technique  and  the  Cross-Multiplication  technique.  With  this 
tool,  the  designer  can  determine  system  stability.  He  may  then  specify 
system  performance  by  pole  placement.  Finally,  he  may  verify  that  performance, 
in  the  time  domain,  with  the  Cross-Multiplication  technique,  for  any  specified 
deterministic  reference  input  the  designer  might  wish  to  select.  This  has 
all  been  programmed  in  a  digital  computer  to  ease  the  computational  burden 
of  the  designer.  Indeed,  the  two  theoretical  methods  have  been  deliberately 
cast  in  a  form  that  is  readily  amenable  to  digital  implementation. 

Although  long  term  in  nature,  the  efforts  of  the  Guidance  and  Control 
Analysis  Group  that  were  supported  by  the  Control  Dynamics  Company  are  going 
along  smoothly.  The  usual  "blind  alleys"  were  encountered  (but  eventually 
circumvented!)  as  the  model  of  the  plant  was  being  defined. 

The  PERSHING  II  G&C  efforts  were  funnel ed  into  a  study  of  the  autopilot 

« 

design  being  implemented.  A  number  of  recommendations  have  been  offered, 
where  they  seemed  appropriate. 

An  effort  to  evaluate  the  consequences  of  microprocessor  implementation 
of  digital  filters  was  begun  in  June  1980.  This  effort  was  interrupted  by 
the  emergence  of  added  emphasis  of  ATBM  G&C. 
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Finally,  a  Digital  Design  Seminar  was  conducted  for  selected  members 
of  the  Guidance  and  Control  Directorate.  It  appears  that  those  whose 
perserverance  permitted  found  the  techniques  they  learned  to  be  useful  and 
the  seminar  to  be  technically  rewarding. 


SECTION  VIII.  RECOMMENDED  FUTURE  EFFORTS 

The  technical  efforts  completed  in  this  contracting  period  have  shed 
light  on  new  technical  work  that  should  be  undertaken.  These  efforts  are 
natural  extensions  of  the  work  that  is  reported  upon  in  this  Final  Report 
and  are  in  three  general  categories: 

A.  Development  of  Discrete  Technology 

In  Section  III  of  this  report,  the  Parameter  Space  technique  was 
expanded.  Several  powerful  future  extensions  to  the  technique  became 
apparent  as  this  work  progressed.  First,  it  probably  can  be  extended  to 
the  popular  w-domain,  but  even  more  important,  the  technique  would  be  even 
more  powerful  if  a  means  of  obviating  the  requirement  for  the  characteristic 
equation  were  developed.  These  two  thrusts  toward  extending  the  technique 
should  be  investigated  and,  if  feasible,  implemented.  A  third  thrust  that 
became  apparent  is  the  need  to  incorporate  modern  sampled-data  techniques 
into  the  present  optimal  G&C  techniques,  thereby  investigating  and  develop¬ 
ing  techniques  for  optimizing  digital  G&C  systems. 

B.  Design  of  Advanced  Discrete  G&C  System 

In  Section  II  of  this  report.  Advanced  G&C  Systems  were  investigated. 
The  major  thrust  was  toward  the  development  of  a  simulation  tool  with  which 
to  test  and  evaluate  various  candidate  G&C  systems  and  their  subsystems.  In 
Section  VI  of  this  report,  advanced  digital  filters  were  investigated 
briefly.  The  results  of  these  two  efforts  focus  toward  the  need  for  the 
design  and  development  of  an  advanced  discrete  (digital)  G&C  system  that  uses 
advanced  digital  filter  techniques.  Such  a  system  would  be  implemented  on 
microprocessors,  where  possible. 

C.  Evaluation  of  PERSHING  II  Digital  G&C  System 

As  the  design  of  the  PERSHING  II  G&C  system  continues,  the  Control 
Dynamics  Company's  unique  digital  design  and  analysis  expertise  and  experi¬ 
ence  would  be  applied  to  evaluate  the  developing  G&C  system.  When  deemed 
appropriate,  modification  to  that  system  would  be  recommended  with  the  aim 
of  obtaining  improved  system  performance  and/or  decreased  system  complexity/ 
cost. 
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SYSTEM  ISSUES 
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4.  Redundancy/autonomy 

5.  Reaction  time  capability 

b.  Battle  management  system 

1.  Saturation  level  &  performance  near  saturation 

2.  Field  worthiness 

3.  Optimum  utilization  of  resources 

4.  Redundancy 

5.  Kill  assessment 
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7.  What  is  "Achilles’  heel"? 

8.  Requirement  for  manual  participation 

9.  Threat  discrimination  capability 

10.  Reaction  time  requirement 

11.  Threat  prioritization 


c.  Reassess  Navy  threat  scenario  for  sensitivity 


II.  GUIOANCE,  NAVIGATION  AND  CONTROL 
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3.  Performance 
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7.  Sensitivity 
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Abstract 

In  1978  the  US  Army  Missile  Command  embarked 
upon  a  task  to  develop  an  advanced  guidance  and  con¬ 
trol  system  for  future  missiles.  It  was  Intended 
to  "leapfrog"  systems  currently  under  development  in 
order  to  meet  the  stringent  demands  and  constraints 
Imposed  by  targets  with  predicted  characteristics  of 
the  1990s  and  beyond  and  by  the  predicated  battle¬ 
field  environment  of  that  time. 

In  this  paper  the  problem  is  redefined  and  the 
latest  development  program  presented.  Results  that 
have  been  achieved  to  date  are  described,  particu¬ 
larly  in  the  areas  of  mathematical  models  of  the 
missiles  and  their  guidance  and  control  (GA C)  sys¬ 
tems  being  used  for  analyses  and  simulations;  aero¬ 
dynamics;  propulsion;  guidance  laws  being  developed 
and  analyzed;  status  of  development  of  Disturbance 
Accommodating  Control;  signal  processing  to  locate 
and  track  the  target(s);  and  digital  design  tool 
development.  The  paper  is  concluded  with  a  section 
of  future  plans,  to  Include  contractor  support. 

I.  Introduction 

The  US  Army  Missile  Command  {MI COM )  recently 
began  a  task  to  develop  an  advanced  GAC  system  for 
Future  Army  Modular  Missiles.  The  intent  is  to 
"leapfrog"  systems  currently  under  development. 

The  purpose  of  this  paper  is  to  describe  the  work 
that  has  been  completed  within  this  new  task  and  to 
provide  an  indication  of  future  efforts  that  are 
now  planned. 

The  first  step  in  implementing  this  task  was 
to  conduct  a  literature  survey  to  establish  a  tech¬ 
nology  base  starting  point.  Following  this  survey, 
guidance  laws  were  placed  in  five  categories  and 
defined  mathematically.  The  implementation  and 
predicted  performance  of  each  category  was  then  in¬ 
vestigated  and  compared  in  light  of  current  and 
predicted  hardware  and  software  capabilities.  ' 

This  work  was  subsequently  updated  in  1979.  2 

The  program  objectives  are  three-fold; 

•  To  develop  and  prove  a  GAC  system  that  is 
capable  of  guiding  and  controlling  future  US  Army 
missiles  (generally  defined  as  air  defense  and  sur¬ 
face-to-surface  general  support)  to  destroy  pre¬ 
scribed  lines  of  future  targets.  This  must  be 
accomplished  under  the  predicted  severe  battlefield 
environment  of  the  future. 

•  To  "leapfrog"  systems  under  current  develop¬ 
ment. 


*Thi s  paper  is  declared  a  work  of  the  US  Govern¬ 
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•  To  broaden  and  deepen  the  existing  GAC  sys¬ 
tem  technology  and  design  base  within  the  GAC 
Directorate  of  MICOM. 

In  the  sequel  the  expected  threat  and  its 
characteristics  will  be  summarized.  This  is  follow¬ 
ed  by  a  description  of  the  development  plan  for  the 
advanced  GAC  system.  The  progress  to  date  is  pre¬ 
sented,  and  the  paper  is  concluded  with  a  section 
on  future  plans. 

II.  Expected  Threat 

Theater  defense  typically  is  provided  by  a 
mixture  of  ground-based  and  airborne  defense  systems 
supported  by  radars,  command  and  control  systems, 
electronic  warfare  equipment,  and  passive  measures 
such  as  camouflage,  decoys,  and  equipment  disper¬ 
sion.  The  air  defense  objective  of  ground  based 
systems  is  to  limit  the  opponent's  effectiveness  by 
attacking  his  critical  assets  so  that  land  forces 
may  maneuver  with  a  minimun  of  interference  from 
the  enemy  air  weaponry. 

For  many  years  now,  the  enemy  doctrine  has 
emphasized  large  mass  and  brute  force,  and  his  air 
attacks  will  provide  no  exception.  It  is  entirely 
feasible  to  assume  that  an  attack  in  the  Central 
Europe  area  will  ne  accompanied  by  several  thousand 
combat  aircraft.  3  in  addition,  his  doctrine  calls 
for  the  massing  of  large  quantities  of  artillery 
fire  on  a  section  selected  for  a  tank-led  break¬ 
through.  It  is  unlikely  that  NATO  forces  either 
now  or  in  the  near  future  will  match  the  Warsaw 
Pact  forces  in  terms  of  numbers  of  weapons,  nor  is 
it  the  intent  to  aim  toward  that  end.  Rather,  it  is 
important  to  optimize  the  effectiveness  of  our, 
smaller  force  to  meet  the  anticipated  threat. 

The  Army  is  attempting  to  maximize  the  effec¬ 
tiveness  of  its  current  family  of  air  defense 
weapons  while  concurrently  developing  a  new  family 
to  meet  the  threat  of  the  1990’s.  In  the  near 
term,  there  will  be  continued  modification  of  cur¬ 
rent  systems  as  necessary,  and  while  still  feasible, 
to  overcome  qualitative  and  quantitative  deficien¬ 
cies.  Longer-term  replacements  continue  in  devel¬ 
opment  or  procurement  for  all  the  major  field  army 
air  defense  systems.  Examples  of  this  strategy 
fnclude  the  following:  high  to  medium  altitude 
missile  systems  -  PATRIOT  for  NIKE  HERCULES  and 
HAWK;  short  r»nge  missile  systms  -  U.  S.  ROLAND  for 
CHAPARRAL;  ^portable  missiles  -  STINGER  for 
REDEYE;  mobile  gun  systems  -  DIVAD  given  to  VULCAN. 
The  systems  will  provide  the  effective  aerial 
umbrella  needed  by  our  forces  to  not  only  survive, 
but  to  fight  effectively. 

For  security  reasons,  it  is  impossible  in  this 
forum  to  describe  specifically  the  air  threat  that 
will  be  encountered  in  the  scenario  described  above. 
However,  in  order  to  quantify  the  problem  somewhat, 
we  shall  attempt  to  attribute  vehicle  characteris¬ 
tics  to  the  enemy  based  on  our  current  technology 
in  the  field  of  air  defense  targets.  The  Targets 
Management  Office,  US  Army  Missile  Command  since 
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1964  has  published  an  extensive  library  of  target 
program  reports.  In  particular,  they  have  classi¬ 
fied  a  variety  of  aerial  targets  as  test  and  evalua¬ 
tion  (T&E)  targets  used  for  air  defense  weapon 
systems.  Needless  to  say,  targets  used  for  this 
purpose  must  exercise  an  air  defense  weapon  system 
to  the  limits  of  its  capabilities.  ' 

A  particularly  Interesting  target  is  known  as 
HAHST,  an  acronym  for  High  Altitude  High  Speed 
Target.  It  is  designed  to  achieve  speeds  up  to 
Mach  4  at  altitudes  up  to  100,000  feet.  Additional 
performance  characteristics  for  HAHST  and  other 
existing  targets  are  given  in  References  2  and  5. 

The  effectiveness  of  any  missile  system  is  con¬ 
ditioned  on  its  ability  to  function  in  an  Electronic 
Countermeasures  environment.  Stand  off  jammers, 
barrage  jammers  and  even  dispensed  chaff  will  be 
used  to  deny  the  air  defense  tracking  radars  the 
capabilities  needed  to  be  effective  in  their  sector. 
The  jammers  are  intended  to  reduce  the  acquisition 
range  of  the  radars  and,  if  perfected,  will  elimi¬ 
nate  accurate  tracking  entirely.  Additionally,  the 
threat,  aircraft  and  missiles  will  be  fabricated  to 
present  the  smallest  possible  radar  cross  section. 
The  state-of-the-art  in  this  field  is  beyond  the 
scope  of  presentation  in  this  paper. 


utilize  different  modules  such  as  propulsion,  guid¬ 
ance,  warhead,  etc. 

A  program  plan  initiated  by  the  G&C  Director¬ 
ate,  MIC0M,  was  undertaken  approximately  two  years 
ago  with  the  program  objectives  enumerated  above. 
The  program  plan  contains  six  intermediate  objec¬ 
tives,  the  end  of  which  each  constitutes  a  program 
milestone.  They  are: 

•  Define  overall  program; 

e  Collect  elements  that  will  form  the  candi¬ 
date  G&C  systems; 

e  Define  candidate  G&C  systems; 

e  Evaluate  the  candidate  G&C  systems: 
select  best  system; 

•  Design  and  fabricate  the  selected  system; 

•  Demonstrate  "Proof-of-Concept"  of  the 
selected  system. 

These  separate  program  elements  (or  intermediate) 
objectives  are  described  below  in  more  detail.  An 
accompanying  milestone  chart  is  provided  as  Fig.  1. 


From  a  defensive  viewpoint,  the  effect  of 
enemy  jamming  of  the  air  defense  radars  and  the 
minimization  of  enemy  attack  aircraft  radar  cross 
sections  have  a  profound  impact  on  the  air  defense 
missile  system.  The  rationale  is  reasonable  and 
straightforward.  If  the  enemy  does  indeed  have 
aircraft  and  attack  missiles  and  remote  pilotless 
vehicles  (RPVs)  either  with  or  better  than  the 
characteristics  attributed  to  them  via  the  method 
above,  and  if  the  enemy  minimizes  his  radar  cross 
section  to  the  current  state-of-the-art,  his  attack 
vehicles  will  be  extremely  difficult  to  acquire  at 
long  range.  The  effect  of  the  combined  high  speed 
and  high  agility  (i.e.,  high  g-maneuver  capability) 
with  low  radar  cross  section  yields  precious  little 
reaction  time  to  the  air  defense  system.  The  close- 
in  acquisition  will  seriously  degrade  the  existing 
air  defense  missile's  G&C  system's  performance, 
since  most  are  based  on  a  proportional  navigation 
and  guidance  (PNG)  law.  An  environment  such  as 
presented  above,  however,  can  be  better  addressed 
in  terms  of  guidance  laws  explicitly  tailored  to 
this  type  of  threat.  Thus,  the  optimally  guided 
and  controlled,  highly  maneuvering,  defensive  mis¬ 
sile  using  terminal  guidance  sensors  chosen  from 
across  a  wide  range  of  the  frequency  spectrum  must 
be  initiated  into  the  development  cycle.  As  a 
necessary  first  step,  the  research  described  in  the 
sequel  addresses  that  problem. 

I 

III.  Development  Program  Plan 

From  an  overall  systems  viewpoint,  th,is  pro¬ 
gram  shall  address  the  issue  of  creating  new  theory 
in  the  G&C  area  to  meet  the  high  performance  threat 
of  the  future  as  one  of  the  leading  technology 
items.  Closely  associated  with  it  and  in  parallel 
with  the  G&C  effort,  weapon  system  work  shall  be 
undertaken  to  modify  airframe  and  propulsion  to  be 
capable  of  engaging  the  threat  of  the  1990's.  Gen¬ 
eral  support  weapons  shall  be  viewed  initially  as  a 
subset  of  the  air  defense  system(s),  whereas  pre¬ 
viously,  these  two  classes  of  weapons  were  develop¬ 
ed  Independently.  This  research  shall  attempt  to 
view  them  as  potentially  similar  systems  that 
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Figure  1  .  Program  Milestone  Chart. 
Task  1.  Oefine  Overall  Program 


In  this  first  task  it  is  necessary  to  define 
the  program  objectives,  goals,  and  constraints. 

This  is  followed  by  a  survey  of  representative  mis¬ 
sile  plants  and  the  characteristics  of  sensors  and 
effectors  that  might  be  used  (either  available  or 
under  development).  Target  dynamics  and  initial 
conditions  must  be  defined,  as  must  missions  for 
the  missile  system(s)  using  the  advanced  G&C  sys¬ 
tems  that  emanates  from  this  program.  Also  includ¬ 
ed  in  the  first  program  element  is  the  beginning  of 
coordination  with  other  missile-developing  agencies 
and  services.  Finally,  the  time  frame  predicted 
usage  c f  the  G&C  system  must  be  defined.  As  des¬ 
cribed  in  this  paper,  most  of  Task  1  has  been  com¬ 
pleted. 

Task  2,  Collect  Elements  for  Candidate  G&C  Systems 

In  this  task  detailed  definitions  of  the 
characteristics  of  candidate  guidance  laws,  candi¬ 
date  autopilots  (control  laws),  state  estimation 
techniques,  state  truncation  techniques,  and  Dis¬ 
turbance  Accommodation  Control  (DAC)  shall  be 
accomplished.  s  Further,  characteristics  of  all 


I 
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expected  disturbances,  both  external  (to  the  mis¬ 
sile)  and  Internal,  must  be  collected  and  evaluated. 
The  need  for  DAC  and  for  state  estimation  trill  be 
evaluated.  Characteristics  of  all  expected  signif¬ 
icant  system  nonl Inearitles  will  be  collected  and 
their  dynamic  Importance  evaluated. 

Task  3.  Oeflne  Candidate  G&C  Systems 

The  components  and  subsystems  that  have  been 
identified  In  the  earlier  tasks  will  now  be  combin¬ 
ed  into  candidate  G&C  systems.  Missile  structural 
dynamics  will  be  determined  (if  not  already  com¬ 
pleted  within  the  plant  definition  of  the  earlier' 
tasks)  or  refined,  as  necessary.  The  simulation 
program  objectives  will  be  defined;  the  simulation 
program  will  then  be  defined;  and  finally,  develop¬ 
ment  of  the  final  system  simulation  will  begin. 
Investigations  of  the  dynamics  of  the  candidate  G&C 
systems  will  begin  using  both  mathematical  analysis 
and  computer  simulation.  A  figure  of  merit  (cost 
functional)  will  be  developed  during  this  task,  as 
will  model  error  effects  and  criteria. 

Task  4.  Evaluate  Candidate  G&C  Systems 

The  development  of  the  system  simulation  and 
the  investigation  of  the  candidate  G&C  systems, 
begun  in  Task  3,  will  be  completed.  The  candidate 
G&C  systems  then  will  be  evaluated  with  respect  to 
the  Figure  of  Merit  developed  in  Task  3,  using 
mathematical  analysis  and  computer  simulation.  The 
strengths  and  weaknesses  of  distributed  versus  cen¬ 
tralized  controllers  will  be  assessed,  probably  as 
enhanced  by  the  use  of  microcomputers .  Finally, 
the  best  G&C  system  will  be  selected  and  integrated 
into  a  missile  airframe. 

Task  5,  Oeslqn  and  Fabrication 

A  detailed  design  of  the  selected  G&C  system 
will  be  performed.  The  testbed(s)  selected  for  use 
in  the  Proof-of-Concept  phase  will  be  fabricated 
and/or  assembled. 

Task  6,  Proof-of-Concept 

The  Proof-of-Concept  will  be  demonstrated  with 
hardware  firings  during  this  task.  These  firings 
will  be  augnented  by  simulations  (computer  and 
hardware-in-the-loop  as  deemed  necessary)  and  anal¬ 
yses  as  necessary. 

IV.  Progress  to  Date 


with  other  organizations  are  anticipated  as  the 
scope  of  the  program  grows . 

The  technical  work  comprising  Tasks  1  and  2 
has  been  apportioned  to  various  members  of  the 
MICOM  and  contractor  team.  As  indicated,  most  of 
Task  1  has  been  completed,  and  effort  on  Task  2 
work  is  underway. 

A.  Target  Definition 

A  comprehensive  investigation  of  predicted 
future  targets  and  their  dynamics  has  been  complet¬ 
ed.  This  Investigation  Included  reviewing  and  dis¬ 
cussing  material  available  within  the  sources  of 
the  US  Army  and  US  Air  Force  and  included  Inputs 
from  several  industrial  organizations.  Finally, 
the  collected  information  was  discussed  with  ele¬ 
ments  of  the  Office  of  the  Under  Secretary  of 
Oefense,  Research  and  Engineering.  While  most  of 
the  collected  information  is  classified  for  secu¬ 
rity  reasons,  the  results  of  this  most  important 
phase  of  the  study  lead  to  emphasizing  targets  in 
three  categories:  highly  maneuverable  aircraft, 
cruise  missiles,  and  tactical  ballistic  missiles. 
Not  addressed  in  this  study  are  RPV's  and  surface 
targets.  It  is  felt  that  use  of  a  sophisticated 
missile  system  to  engage  numerous  RPV's  would  not 
be  cost  effective  and  that  missiles  currently 
under  development  will  be  able  to  combat  surface 
targets.  Except  for  periodic  updating,  this  sub¬ 
task  is  now  complete. 

B.  Missile  Plant  Characteristics 

In  order  to  get  the  other  subtasks  underway, 
standard  equations  of  motion  for  missile  bodies 
with  some  structural  flexibility  have  been  used. 

The  missiles  are  assumed  to  be  acted  upon  by  the 
usual  aerodynamic  forces  and  torques.  Currently, 
variations  of  the  SPRINT  missile  family  are  among 
the  principal  contenders  for  the  airframe  and  pro¬ 
pulsion  system. 

C.  Aerodynamics 

A  detailed  plan  is  being  generated  for  future 
actions  in  gathering  and  generating  aerodynamic 
data  to  be  used  in  this  program.  Novel  aerodynamic 
shapes  are  under  consideration  and  evaluation  for 
use  In  developing  control  authority  for  missiles 
using  advanced  G&C  systems.  These  include  missiles 
using  configurations  shown  in  Table  1  (along  with 
their  characteristics). 


The  accomplishments  to  date  have  been  achieved 
by  elements  of  MICOM  and  through  the  use  of  re¬ 
search  contracts.  It  has  been  supplemented  through 
coordination  with  the  US  Air  Force  Armaments  Labo¬ 
ratory  at  Eglln  Air  Force  Base,  the  US  Army  Bal¬ 
listic  Missile  Oefense  Systems  Command  at  Hunts¬ 
ville,  Alabama,  and  the  Office  of  the  Under  Secre¬ 
tary  of  Oefense.  This  coordination  has  been  and 
will  continue  to  be  carried  out  to  eliminate  dupli¬ 
cated  development  effort  on  similar  projects  within 
the  Department  of  Defense.  To  augnent  the  in-house 
research  and  engineering  capability,  several  re¬ 
search  contracts  have  been  initiated  in  specialized 
areas  currently  including:  the  Computer  Sciences 
Corporation,  Huntsville,  Alabama;  the  Oynamic 
Systems  Research  and  Training  Corporation,  Hunts¬ 
ville,  Alabama;  the  University  of  Florida;  Western 
Kentucky  University;  and  the  Control  Dynamics  Com¬ 
pany,  Huntsville,  Alabama.  Additional  contracts 


0.  Error  Sources 

Expected  error  sources  have  been  categorized 
into  five  detailed  groups,  to  include  predicted 
mean  and  standard  deviation  values.  Several  are 
presented  as  representative  of  the  groups.  These 
will,  of  course,  change  as  the  choices  are  narrowed 
and  the  program  develops. 

1.  Missile.  Included  in  this  group  are  para¬ 
meter  uncertainties  such  as  thrust  magnitude  and 
misalignment,  C.G.  location  and  offset,  mass,  and 
transverse  and  axial  moments  of  inertia. 

2,  Aerodynamic.  Normal  and  axial  forces; 
pitch,  yaw,  and  roll  moments;  pitch  and  yaw  damping 
derivatives;  roll  damping  coefficient  due  to  fin 
deflection;  and  axial  drag  are  factors. 


Table  1.  Aerodynamic  Configurations  and  Characteristics. 


Configuration 

Advantages 

Disadvantages 

Mach 

Mo. 

Altitude 

(ft.) 

Maneuverability 

(g‘s) 

1 .  Deployable  wings 
with  all-moveable 
tails 

a.  Simple  con¬ 
trols 

b.  Good  roll 
control 

a.  Control  force 

In  opposite  direc¬ 
tion  from  maneuver 

b.  High  angle  of 
attack 

2-4 

0-100  K 

10-15 

2.  Low  aspect 
ratio,  long-chord 
delta  wing  with 
tab  controls 

a.  Simple  con¬ 
trols 

b.  Good  roll 
control 

a.  Control  force 

In  opposite  direc¬ 
tion  from  maneuver 

b.  High  angle  of 
attack 

4-7 

0-140  K 

20-30 

3.  Flared  skirt- 
stabilized  missile 
with  all  moveable 
wings 

a.  Control  force 
in  direction  of 
maneuver 

b.  Lower  angle 
of  attack 

Possible  higher 
drag  and  hinge 
moment 

• 

4-7 

0-140  K 

-i 

20-30 

4.  Reentry  body 
shape  wltn  moveable 
wedge  controls  * 

a.  Simple  shape 

b.  Good  roll 
control 

Roll  control  sur¬ 
faces  separate 
from  maneuvering 
controls 

5-10 

7 

? 

5.  Lifting  body  * 

a.  Lower  angle 
of  attack 

b.  Control  force 
In  direction  of 
maneuver 

Complicated  aero¬ 
dynamic  shape 

3-6 

0-90  K 
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the  missile;  target  velocity  and  illumination 
jitter;  semiactive  laser  pointing  accuracy  and  beam 
divergence;  gravity  bias;  atmospheric  effects 
(e.g.,  upon  radio  range). 

5.  Subsystems 

•  Common  nonlinearities  -  saturation, 
coulomb  (and  other)  friction,  backlash,  and  bang- 
bang  with  dead  zone. 

•  Computer  -  quantization,  truncation,  and 
fixed  word  length. 

•  Seeker  -  boresight  error  (in  pitch  and 
yaw)  due  to  servo  noise;  channel  crosscoupling; 
coupling  between  the  seeker  head  and  the  airframe; 
effects  of  the  radome  and  irdome  on  angle  linear¬ 
ity;  angle  bias  (the  electrical  equivalent  of  mech¬ 
anical  BSE);  gain  stability;  angle  noise;  and  bore- 
sight  error  in  pitch  and  yaw  due  to  clutter,  re¬ 
ceiver,  and  jamming  noise. 

e  Autopilot  -  bias  errors,  time  delays, 
and  gain  stability. 

•  Guidance  -  errors  in  Initial  position, 
velocity,  and  acceleration. 

As  the  program  progresses,  these  error  sources 
will  be  analyzed  further  to  determine  which  might 
be  amenable  to  cancellation  by  appropriate  disturb¬ 
ance  accomodation  design  theory  and  techniques. 


*  Conenent.  Must  use  bank-to-turn  guidance. 


3.  Instruments 

•  Accelerometers  -  scale  factor  stability, 
bias  stability,  non-orthogonality,  "g2"  -  scale 
factor,  third-order  scale  factor,  cross-axis  sensi¬ 
tivity,  cross-coupling,  scale  factor  asymmetry,  and 
rectification  error. 

•  Gyroscopes  -  scale  factor  stability; 
bias  stability;  non-orthogonality;  anisoelastic 
drift;  drift  rates  in  pitch,  yaw,  and  roll  due  to 
input  axis  and  spin  axis  accelerations  as  well  as 
those  (rates)  Independent  of  acceleration,  due  to 
torquer  nonl inearlties,  and  due  to  electronic 
noise;  and  mass  unbalance. 

•  Porro-prism  azimuth  alignment. 

e  Laser  inertial  measurement  unit  mis¬ 
alignment  with  respect  to  missile  body  axes. 

e  Accelerometer  triad  origin  displacement. 

•  Uplink/downlink  bias  and  calibration 


errors . 


•  Optical  correlator  errors. 


e  Radar  -  range  track  and  angle  track 
noise  and  accuracy,  ground  clutter  noise,  and  tar¬ 
get  glint  noise. 

4 ,  External  ■  Wind  magnitude  and  direction 
(initial  azimuth,  elevation,  and  roll  alignment  of 


t.  Guidance  Laws 

This  Is  a  major  thrust  area  within  the  program. 
As  indicated  above,  an  extensive  literature  search 
has  been  completed  and  documented.  This  was  follow¬ 
ed  by  placing  guidance  laws  in  five  categories  and 
describing  each  mathematically.  The  implementation 
and  predicted  performance  of  each  category  has  been 
initially  reported  in  Reference  1  and  subsequently 
updated  in  Reference  2.  A  sunmary  of  the  latter 
is  shown  in  Table  2. 

Investigations  continue  into  areas  of  optimal 
guidance,  in  particular,  emphasizing  digital  as¬ 
pects  in  anticipation  of  the  expected  use  of  on¬ 
board  digital  controllers.  Constant  effort  is  made 
to  reduce  implementation  complexity,  where  complex¬ 
ity  is  defined  as  requirements  for  hardware  and 
software.  Other  innovative  techniques  are  under 
investigation.  They  include  assessing  the  potential 
application  of  Singular  Perturbation  Theory,  Distur¬ 
bance  Accommodating  Theory,  and  means  of  determining 
or  predicting  the  very  important  (for  optimal 


applications)  quantity,  time-to-go,  i.e.,  remaining 
time  of  flight  at  any  instant  (see  Table  3). 

1 .  Optimal  Control .  A  conventional  implementa 
tfon  of  Linear-Quadratic  Optimal  Control  is  por¬ 
trayed  graphically  in  Figure  2.  The  optimal  con¬ 
trol  authority,  u®q,  is  selected  to  be 

u°Q  -  -R*Vkxx  .  (1) 

The  matric  Riccati  equation  is  solved  to  obtain  the 
control  gain,  <x ( t ) ,  while  minimizing  a  quadratic 

performance  index. 

2,  Disturbance  Accomodating  Control .  The 
Disturbance  Accommodating  Control  (OAC)  theory  is 
described  for  the  continuous-time  domain  in  Refer¬ 
ence  6.  The  general  nature  of  the  DA C  controller 
is  to  generate  a  real-time  on-line  estimate  of  the 
actual  (instantaneous)  disturbance  waveform  and 
create  a  special  control  action  that  exactly 


Table  2.  Conventional  Guidance 


Approach 

Advantages 

Limitations 

1.  Attitude  Pursuit 

a.  Simplest  implementation 

b.  Fixed  targets 

Sensitive  to  target  velocity, 
disturbances 

2.  Velocity  Pursuit 

a.  Simple  implementation 

b.  Non-maneuvering  targets 

Sensitive  to  target  accelera¬ 
tion,  disturbances 

3.  Proportional  Naviga¬ 
tion  Guidance  (PNG) 

a.  Simple  implementation 

b.  Maneuvering  targets 

Sensitive  to  high  end  game 
maneuvers 

Table  3.  New  Guidance  Approaches 


Approach 

Advantages 

Limitations 

1.  Linear  Quadratic  (LQ) 
Regulator 

Better  than  PNG  against 
maneuvering  target 

a.  Tqq  estimate  required 

b.  Must  compute  time-varying 
control  gain  (1^) 

c.  Disturbances  ignored 

2.  Linear  Quadratic 

Gaussian  (LQG) 

Regulator 

Better  than  LQ  against 
noise* type  disturbances 

a.  Tgg  estimate  required 

b.  Must  compute  time-varying 
control  gain  (1^) 

c.  Must  compute  Kalman  gain 
for  estimator 

3.  Disturbance-Utiliz¬ 
ing  Control  (DUC) 

Better  than  LQ  or  LQG 
against  waveform- type 
disturbances 

a.  Tgg  estimate  required 

b.  Must  compute  2  time-vary¬ 
ing  control  gains  (Kx,  K^) 

4.  Singular  Perturbations 

Computational  efficiency 

An  approximation  to  optimal 
control 

cancels-out  the  disturbance  effect  on  the  missile. 
The  DAC  theory  will  be  extended  Into  a  discrete¬ 
time  domain.  To  date,  the  class  of  systems  and 
disturbances  amenable  to  a  discrete-time  version  of 
OAC  have  been  defined.  Work  Is  underway  on  a  des¬ 
cription  of  the  state-reconstructor  that  will  be 
associated  with  a  digital  DAC  controller.  Also 
under  consideration  is  the  possibility  of  applying 
OAC  to  various  missile  subsystem  or  component  out¬ 
puts,  such  as  sensors,  whose  normal  outputs  have 
been  modified  by  the  Influence  of  the  disturbances 
upon  the  sensors. 

An  Innovative  modification  of  DAC  is  Dis¬ 
turbance-Util  izlng  Optimal  Control  (OUC)7.  In  this 
case  waveform-type  disturbances  are  exploited  opti¬ 
mally.  Examples  of  such  disturbances  are  drag, 
target  maneuvers,  wind  gusts,  any  effects  of  the 
gravitational  field.  A  graphical  portrayal  of  DUC 
implementation  Is  shown  In  Figure  3.  In  this  case, 
the  optimal  control  authority  Is  specified  as 

u°0UC  *  -R'V(Kx*  +  Kxz*>  .  (2) 

where  *x(t)  is  found  by  solving  a  matric  Riccati 
equation,  and  K  (t)  Is  found  by  solving  a  linear 
differential  equation.  To  date,  computer  simula¬ 
tions  have  shown  OUC  to  be  quite  effective  when 
compared  to  the  performance  achieved  by  using  con¬ 
ventional  IQ  controllers. 


3.  Singular  Perturbations .  Appl  icatlon  of 
Singular  Perturbation  Theory  to  missile  control  may 
be  attractive  if  It  Is  deemed  necessary  for  the  g 
control  law  to  account  for  high  order  model  terms. 
The  standard  approach  is  to  approximate  the  model 
with  relatively  lower  order  equations.  However, 
the  neglected  higher  order  terms  may  be  dynamically 
significant.  Because  their  Inclusion  might  create 
computational  problems,  a  possible  alternative 
might  be  the  application  of  Singular  Perturbation 
Theory.  The  application  is  particularly  amenable 
to  controller  design  where  the  open-loop  plant  has 
a  wide  eigenvalue  dispersion,  slow  and  fast  modes, 
or  parasitic  parameters. 

F,  Sensor  Characteristics 

Various  existing  and  predicted  sensors  have 
been  characterized.  Trade  studies  to  aid  in  their 
selection  have  been  identified.  It  now  appears 
that  a  sensor  using  some  form  of  pattern  recogni¬ 
tion  may  be  required.  The  problem  is  to  find  or 
begin  developing  a  sensor  that  can  provide  a  guid¬ 
ance  signal  with  a  superior  signal -to-nolse  ratio 
from  data  that  has  been  deliberately  modified,  for 
example,  by  high  powered  jamming  equipment. 

G.  Computer-Aided  Design  Tools 

A  number  of  computer  programs  are  available 
(such  as  root  locus  and  other  frequency  domain 
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Figure  2.  Conventional  Linear  Quadratic  Optimal  Control. 


Disturbance 


Figure  3.  Disturbance-Utilizing  Optimal  Control. 
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techniques)  to  aid  In  design  of  control  systems. 
Most  of  them  are  In  the  continuous-time  domain, 
although  there  does  exist  a  z-transform  manipula¬ 
tion  program.  An  executive  routine  currently  Is 
being  coded  to  manipulate  these  programs  effi¬ 
ciently. 

H.  Microcomputer  State-of-the-Art 

A  continuing  assessment  of  microcomputer  state- 
of-the-art  has  been  Instigated  within  this  program. 
Microcomputers  will  become  an  integral  part  of 
advanced  G&C  systems,  with  their  small  size,  light 
weight,  and  relatively  low  cost.  Presently,  micro¬ 
computers  software  capabilities  are  being  Investi¬ 
gated  with  a  Tektronix  8002  Microprocessor  Labo¬ 
ratory.  This  Investigation  will  be  followed  by  an 
analysis  of  the  computational  capability  to  meet 
the  requirements  imposed  by  an  advanced  G&C  system. 

I.  Digital  Design  Tool  Development 

The  probability  that  an  advanced  G&C  system 
will  be  implemented  digitally  seems  to  be  nearly 
assured.  It  is  for  this  reason  that  the  foregoing 
work  has  been  oriented  strongly  in  the  direction 
of  digital  implementation  and  its  consequences.  A 
review  of  analysis  and  design  tools  available  to 
the  G&C  system  designer  is  being  conducted.  Where 
the  need  for  a  new  tool  appears  to  be  warranted,  it 
will  be  developed  within  this  program  when  possible. 
Emphasis  is  placed  on  both  simplicity  of  applica¬ 
tion  and  being  able  to  draw  on  practising  engineers 
existing  engineering  training  and  experience. 

One  example  of  such  a  tool  is  "SAM"  (acronym 
for  Systematic  Analytical  Method).  SAM  provides  an 
alternative  to  the  use  of  signal  flow  graphs  and 
the  application  of  Mason's  Gain  Rule  to  determine 
selected  states  of  the  missile.9  The  technique  is 
particularly  useful  for  analysis  of  complicated 
sampled -data  control  systems.  An  advantage  of 
using  SAM  is  that  the  cumbersome  application  of 
Mason's  Gain  Formula  can  be  avoided.  Further,  the 
entire  method  of  constructing  signal  flow  graphs 
may  be  circumvented.  Since  only  the  equations  des¬ 
cribing  the  system  are  needed  for  SAM,  even  the 
customary  block  diagram  is  not  needed.  The  tech¬ 
nique  is  analytical  in  nature  and  makes  use  of  a 
systematic  manipulation  of  the  system  algebraic 
equations.  These  manipulations  follow  prescribed 
rules  set  forth  in  the  technique. 

A  second  simplified  technique  is  the  determina¬ 
tion  of  digital  control  system  response  by  cross- 
multiplication.10  This  technique  permits  the 
analyst  to  obtain  the  response  (at  the  sampling  in¬ 
stants)  of  any  system  state  from  its  closed-loop 
transfer  function  expressed  in  the  complex  z-domain. 
If  it  is  desired  to  know  the  response  between 
sampling  instants,  either  the  submultiple  method  or 
the  modified  z-transform  method  may  be  adapted  to 
the  cross-multiplication  technique. 

A  third,  more  sophisticated  technique  is  the 
Parameter  Space  Method.  It  is  being  developed  for 
determining  the  stability  and  dynamic  characteris¬ 
tics  of  a  digital  control  system  in  terms  of  sever¬ 
al  selected  system  parameters."  The  method  re¬ 
quires  that  the  system  characteristic  eouatfon  be 
available  in  the  complex  z-domain.  Although  not 
necessary,  its  application  is  facilitated  by  aug¬ 
menting  the  analytical  results  with  graphical  por¬ 
trayals  in  a  selected  multiparameter  space.  The 


method  is  based  on  analysis  and  synthesis  methods 
for  linear  and  nonlinear  control  system  design 
which  are  amply  described  in  Siljak's  excellent 
monograph  on  the  subject. '2  In  essence,  the  para¬ 
meter  space  method  permits  the  designer  to  evaluate 
graphically  the  effects  of  the  locations  of  the 
roots  of  the  characteristic  equation.  Hence,  he 
may  design  the  control  system  in  terms  of  his  se¬ 
lected  performance  criteria;  e.g.,  absolute  stabil¬ 
ity,  damping  ratio,  settling  time.  He  is  able  to 
see  the  effects  on  the  characteristic  equation 
roots  (and  hence  on  system  dynamics)  of  changing 
several  adjustable  parameters.  The  method  has  been 
extended  to  portray  the  effect  of  varying  the 
sampling  period,  thereby  permitting  one  to  observe 
the  effect  of  the  choice  of  various  values  assigned 
to  the  sampling  period  on  absolute  and  relative 
stability.  Also,  simple  recursive  formulas  have 
been  derived  so  that  the  resulting  formulation  is 
deliberately  cast  in  a  form  particularly  amenable 
to  solution  by  a  digital  computer  or  a  desk  cal¬ 
culator,  emphasizing  the  interplay  between  analysis 
and  computing  machines. 

J.  Documentation 

A  major  portion  of  the  progress  reported  above 
has  been  documented  in  a  comprehensive  US  Army 
Missile  Command  report. 13 

V.  Future  Plans 

Plans  for  the  next  fiscal  year  revolve  primar¬ 
ily  about  impl ementating  Task  3  (define  candidate 
G&C  systems)  of  the  development  program  as  well  as 
completing  any  partially-completed  portions  of  Task 
2  now  underway.  Close  communication  will  be  main¬ 
tained  with  the  intelligence  community  to  become 
aware  of  any  changes  to  the  presently  predicted 
targets  as  contrasted  to  the  attribution  mentioned 
above.  Detailed  analytical  models  of  the  missile 
plant(s),  effectors,  and  autopilot(s)  will  be  devel¬ 
oped.  The  possibility  of  modifying  plant  charac¬ 
teristics  by  making  innovative  use  of  aerodynamics 
is  to  be  investigated.  Projected  advances  in  the 
field  of  propulsion  will  also  be  investigated.  The 
comparison  of  various  digitally-implemented  guid¬ 
ance  laws  will  be  continued.  Including  those  that 
incorporate  DAC.  Applications  (other  than  to  guid¬ 
ance  laws)  of  DAC  theory  to  improve  system  perform¬ 
ance  will  be  investigated.  Development  of  a  mod¬ 
ular  guidance  simulation  to  implement  this  investi¬ 
gation  has  already  begun.  Trade  studies  concerning 
identified  existing  and  future  sensors  will  be  con¬ 
ducted,  with  the  possibility  of  developing  a  new 
sensor  with  capabilities  not  yet  in  existence.  Ad¬ 
vances  in  microcomputer  state-of-the-art  will  be 
watched  closely.  Finally,  tools  to  aid  the  G&C 
system  designer  to  handle  digital  implementation 
will  continue  to  be  both  assessed  and  developed. 

VI.  Conclusions 

As  indicated  within  this  paper,  there  is  a  clear 
need  for  the  development  of  an  advanced  G&C  system 
for  US  Army  future  tactical  missiles.  This  need  is 
dictated  by  the  predicted  targets  of  the  future, 
the  anticipated  battlefield  of  the  future,  and  the 
characteristics  of  tactical  missiles  that  are  ei¬ 
ther  in  the  inventory  or  under  development  at  this 
time.  The  nature  of  the  future  G&C  system  will  be 
digital  so  that  the  overall  missile  system  may  be 
availed  of  present  and  predicted  advantages  that 
are  implicit  in  digital  controllers.  The  form  of 
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the  guidance  law  probably  will  be  optimal,  since 
the  performance  criteria  that  must  be  minimized 
must  take  more  into  account  than  miss  distance. 
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1.  INTRODUCTION 


The  parameter  space  method  provides  an  analytical  tool  developed  for  use  in  control  system 
analysis  and  synthesis.  Although  not  necessary,  its  application  is  facilitated  by  augmenting  the 
analytical  results  with  graphical  portrayals  in  a  selected  multi-parameter  space.  The  method 
requires  that  the  control  system  be  described  by  a  characteristic  equation  which,  for  sampled- 
data  or  digital  systems,  may  be  expressed  in  the  z-domain.  The  technique  is  based  on  the 
analysis  and  synthesis  methods  for  linear  and  nonlinear  control  system  design  which  are 
amply  described  in  Siljak’s  excellent  monograph  on  the  subject.1  In  another  work,  Siljak 
describes  the  application  of  the  technique  to  the  analysis  and  synthesis  of  linear  sampled-data 
control  systems.2 

Once  the  system  characteristic  equation  has  been  obtained,  the  parameter  plane  method 
enables  the  designer  to  evaluate  graphically  the  locations  of  roots  of  the  equation.  Hence,  he 
may  design  the  control  system  in  terms  of  the  chosen  performance  criteria;  e.g.  absolute 
stability,  damping  ratio,  and  setting  time.  He  is  able  to  see  the  effect  on  the  characteristic 
equation  roots  of  changing  two  adjustable  parameters.  Siljak  further  simplified  the  design  by 
introducing  Chebyshev  functions  into  the  equations,  thereby  putting  them  in  a  form 
that  is  amenable  to  their  solution  by  a  digital  computer. 

The  method  has  been  extended  to  portray  the  effect  of  varying  the  sampling  period.3’4  The 
extended  method  permits  one  to  see  the  effect  of  the  choice  of  values  assigned  to  the  sampling 
period  on  absolute  and  relative  stability.  Also,  the  recursive  formulas  shown  therein  are 
simpler  in  form  than  the  Chebyshev  functions  of  Reference  2.  The  resulting  formulation  is 
deliberately  cast  in  a  form  that  makes  it  particularly  amenable  to  solution  by  a  digital 
computer  or  a  desk  calculator,  again  emphasizing  the  interplay  between  analysis  and 
computing  machines.  When  portrayed  graphically,  the  results  show  the  dynamic  relation 
between  the  selected  parameters  and  the  characteristic  equation  roots,  as  a  function  of  the 
nondimensional  independent  argument,  a»nT.  Hence  one  readily  can  deduce  the  dynamic 
effect  upon  the  system  of  various  combinations  of  values  of  the  selected  parameters  defining 
the  parameter  space. 


1.  D.D.  Siljak,  Sonlinear  Systems.  Wiley,  New  York.  1969. 

2.  D.D.  Siljak,  Analysis  and  Synthesis  of  Feedback  Control  Systems  in  the  Parameter  Plane.  Part  U-Sampled  Dau  Systems." 
Transactions  of  the  Institute  of  Electrical  and  Electronics  Engineers.  Part  II  Applications  and  Industry.  Voi.  83.  November 
1964.  pp.  458-466. 

3.  S.M.  Seltzer.  “Sampled-Data  Control  System  Design  in  the  Parameter  Plane."  Proceedings  of  the  8th  Annual  Allerton 
Conference  Circuits  and  System  Theory.  1970,  pp.  454-463. 

4.  S.M.  Seltzer.  "Enhancing  Simulation  Efficiency  with  AnalyttcalTools."Compurersomf£/errrica/ Engineering,  Vol. 2, 1975, 
pp.  35-44. 
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The  history  of  the  continuous-time  domain  version  of  the  parameter  plane  technique  is  well- 
described  with  suitable  references  in  Reference  1.  Briefly  summarizing  that  history,  in  1876 
I.  A.  Vishnegradsky  of  the  Leningrad  School  of  Theoretical  and  Applied  Mechanics  developed 
and  used  the  first  version  of  the  parameter  plane  technique  to  portray  system  stability  and 
transient  characteristics  of  a  third  order  system  on  a  two-parameter  plane.  In  1949  Professor 
Yu  1.  Neimark  of  the  Russian  School  of  Automatic  Control  generalized  Vishnegradsky's 
approach  to  permit  the  decomposition  of  a  two-parameter  domain  (D)  describing  an  nth  order 
system  into  stable  and  unstable  regions.  The  technique  was  called  D-decomposition.  During 
the  period  1959-1966  Professor  D.  Mitrovic,  founder  of  a  Belgrade  group  of  automatic 
control,  extended  the  method  to  enable  the  analyst  to  relate  the  system's  variable  parameters 
to  the  system  response,  using  the  last  two  coefficients  of  an  nth  order  characteristic  equation. 
Beginning  in  1964.  Professor  D.D.  Siijak,  then  a  student  of  Mitrovic *s  at  the  University  of 
Belgrade,  generalized  the  method  and  called  it  the  Parameter  Plane  method.  His  method 
permitted  the  analyst  to  select  an  arbitrary  pair  of  characteristic  equation  coefficients  (or 
parameters  appearing  within  the  coefficients)  and  portray  both  graphically  and  analytically 
the  dependence  of  the  system  response  upon  the  selected  parameters.  The  method  was 
extended  subsequently  by  Siijak  and  others  to  encompass  a  host  of  related  problems.  In  1967 
Professor  J.  George  modified  the  D-decomposition  method  to  enable  the  portrayal  of  the 
absolute  stability  region  in  a  multi-parameter  space.  George  also  showed  how  to  portray 
contours  of  relative  stability,  as  did  Siijak.  All  of  the  foregoing  work  is  carefully  and 
completely  referenced  within  Reference  1 .  In  1966  and  subsequent  years.  Seltzer  has  applied 
the  parameter  space  method  to  the  design  of  missile,  aircraft,  and  satellite  controllers, 
including  systems  containing  one  or  two  nonlinearities;  the  analysis  of  the  dynamic  effects  of 
the  nonliner  “Solid  Friction”  (Dahl)  model  for  systems  with  ball  bearings,  such  as  control 
moment  gyroscopes  and  reaction  wheels;  and  the  specification  by  the  system  designer  of  the 
dynamic  structural  flexibility  constraints  to  the  structural  designer.  Most  of  the  work  has 
appeared  in  the  technical  journals  of  the  Institute  of  Electrical  and  Electronic  Engineers 
(IEEE),  the  American  Institute  of  Astronautics  and  Aeronautics  (AIAA),  the  International 
Journal  of  Control,  and  the  journal,  Computers  &  Electrical  Engineering.  A  portion  of  the 
history  that  has  not  been  reported  upon  previously  (with  one  exception  to  be  noted)  is  the 
control  system  work  conducted  by  the  German  rocket  scientists  in  the  early  1940’s  at 
Peenemunde.  There,  Dr.  W.  Haeussermann  and  others  applied  the  D-decomposition 
technique  to  the  design  of  the  V-2  Rocket,  following  Dr.  Haeussermann’s  earlier  (pre-World 
War  II)  application  to  the  control  of  an  underwater  torpedo.  This  work  was  not  published  in 
the  open  literature  because  of  national  security  constraints.  When  the  group  came  to  the 
United  States  to  work  with  the  Army  Ballistic  Missile  Agency  (first  at  Ft.  Bliss,  Texas,  then  at 
Redstone  Arsenal,  Huntsville,  Alabama).  Dr.  Haeussermann  and  his  associates  continued  to 
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apply  the  method  to  US  Army  missiles  (and  later  to  NASA  space  vehicles).  Again,  national 
(this  time,  another  nation!)  security  precluded  publication  in  the  open  literature  until  1957.5 

In  1964  Professor  Siljak  published  the  first  application  of  the  parameter  plane  technique  to 
sampled-data  systems.3  As  mentioned  above,  this  was  extended  in  References  3  and  4.  In  197 1 
Seltzer  presented  an  algorithm  for  systematically  solving  the  Popov  Criterion  applied  to 
sampled-data  systems.*  Applications  of  these  sampled-data  parameter  space  techniques  are 
found  in  References  3,  4  and  6. 

2.  ANALYTICAL  PRELIMINARIES 

The  technique  requires  that  the  control  system  be  described  by  a  characteristic  equation  in 
the  z-domain.  Two  adjustable  parameters  (ko.  ki)  are  selected,  and  the  characteristic  equation 
(CE)  is  recast  in  terms  of  them;  i.e. 

n  4 

ce  -  z  yt  z2  -  o,  (i) 

j-0  1 


Yj  “  Yj  (dj*  V  ki} 


Ts  10 
z  *  e  ■  re  , 


(2) 

(3) 


r  -  r(5,wn,T)  -  e'?UJnT’ 

B«  B(£,<jOn>T)  ■  cos0»  cos  »  coswT 


(4) 

(5) 


where  a*  and  T  represent  the  damping  ratio,  natural  frequency,  and  sampling  period, 
respectively.  The  symbol  d,  represents  all  system  parameters  other  than  k«,  ki.  To  transform 
the  characteristic  equation  from  an  nth  order  polynomial  into  an  algebraic  equation,  z!  may  be 
defined  in  terms  of  its  real  and  imaginary  parts,  R,  and  I,: 

zj  -  +  i  Ij  .  (6) 


5.  W  Haeussermann.  “Stability  Areas  of  Missile  Control  Systems."  Jet  Propulsion.  July  I9S7.  pp.  787-795. 

6.  S.M.  Seltzer.  “An  Algorithm  for  Solving  the  Popov  Criterion  Applied  to  Sampled-Data  Systems."  Proceedings  of  the  Third 
Southeastern  Symposium  on  System  Theory.  5-6  April  1971.  pp.  G5-0  through  G5-7. 
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It  readily  follows  from  Equation  (3)  that  the  values  of  R,  and  I,  are 
Rj  *  cos  jO 


(7R) 


and 


1^  -  rj  sin  jQ  .  (71) 

While  equations  (7R)  and  (71)  are  satisfactory  for  determining  the  values  of  Rj  and  lj,  a  set  of 
recursive  formulas  can  be  derived  that  are  particularly  amenable  to  implementation  on  a 
desktop  calculator  or  digital  computer.  They  are  shown  below  as  Equation  (8)  and  are  derived 
in  the  appendix. 

v  -  -  'Vi  ,  «> 

where  Xj  may  be  either  R,  or  Ij.  Only  two  values  of  X,  (i.e.,  two  values  each  of  Rj  and  Ij)  are 
needed  to  obtain  iteratively,  one  at  a  time,  all  other  values  of  R„  I,.  These  are  obtained  from 
the  definition  of  z  in  Equation  (3).  When  j=0,  the  value  of  zj  is  unity  and.  from  Equation  (6), 

z°  -  *0  -  1  xo  <’> 

or 

1  -  R0  (10R) 

0  -  IQ  .  (101) 

When  j—  I,  the  value  of  V  is,  in  Euler  form, 

z*  *  z  *  r  cos  0  +ir  sin  0 

■ynr  (ii) 


Hence, 


and 


Ix  -  rv/l-B2  .  (121) 

It  will  turn  out  to  be  useful  in  the  sequel  to  observe  from  Equation  (8)  that  the  radical,  Vl-B2, 
appears  as  a  factor  in  each  value  of  1,.  Solely  for  simplicity,  it  may  be  factored  out  by  defining  a 
new  term,  i/as 

Ij'  -  Ij/x/l-B2  .  (13) 

It  may  be  observed  that  R,  and  I,  are  functions  of  r  and  B  or  of  £,«>„,  and T.  If  Equation (6)  is 
substituted  into  Equation  ( I )  and  the  real  and  imaginary  parts  of  the  resulting  equation  are 
separated,  two  simultaneous  algebraic  equations  are  obtained. 

These  two  equations  contain  the  adjustable  parameters,  or  variables,  ko,  ki .  Hence  it  may  be 
observed  that  the  two  equations  may  be  solved  explicitly  for  k0  and  ki  as  functions  of  the 
other  system  parameters  (d,)  and,  in  particular,  as  functions  of  the  independent  argument, 
cu„T.  It  is  the  latter  observation  that  forms  the  basis  of  the  parameter  space  method  for 
determining  stability  regions  and  system  dynamic  response  in  terms  of  selected  system 
parameters. 

3.  STABILITY  DETERMINATION 


The  method  involves  the  definition  of  stability  boundaries  on  a  multi-parameter  (to  include 
ko.  ki)  space.  These  boundaries  are  found  from  the  pair  of  simultaneous  algebraic  equations 
that  result  from  the  following  operations  on  the  system  characteristic  equation  written  in 
terms  of  the  complex  variable,  z.  Equation  (6)  is  substituted  into  characteristic  Equation  (1), 
and  the  resulting  real  and  imaginary  parts  are  separated  and  equated  to  zero.  The  two  resulting 
algebraic  equations  may  be  solved  for  any  two  parameters  (such  as  ko  and  ki)  within  them. 

There  may  be  as  many  as  four  stability  boundaries  for  any  system  described  by 
characteristic  Equation  (1),  although  it  is  not  necessary  for  all  four  to  exist. 

One  stability  boundary  separates  the  stable  complex  conjugate  pairs  of  roots  from  the 
unstable  ones.  It  consists  of  a  map  of  the  unit  circle  from  the  complex  z-plane  onto  a  selected 
(such  as  ko.  ki)  parameter  space.  Initially  one  may  consider  the  oft-occurring  case  where  the 
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coefficients  y,  of  the  powers  of  z  in  the  characteristic  equation  are  linear  combinations  of  ko, 
ki,  i.e.. 


Yj  *  aj  ko  +  bj  kl  +  cj  '  (14) 

where  a,,  bj,  c,  represent  all  system  parameters  other  than  ko,  kt.  In  this  case  the  two 
simultaneous  equations  resulting  from  the  real  and  imaginary,  respectively,  parts  of  the 
characteristic  equation  assume  the  form. 

Re  {C.E.}  ■  Aj  s  +  b-C1  -  0,  (15R) 

Im  |C.E. |  -  A2  a  +  B2  b-C2  -  0,  (151) 


and  may  be  solved  for  ko  and  ki  readily  by  applying  Cramer's  Rule,  yielding 


k 

o 


B1C2  ~ 


Vl 


(16) 


k 


1 


A2C1 


(17) 


J  *  A1B2  -  Vl 


(18) 


A  -  E  a  R  Bx  -  E  b  R  Cj  -  E  c  R 

j«o  3  1  j-0  J  \  1  j-0  J 


A2  "  E  a.I,  B2  ■  E  b  I,  C2  ■  E  c.I. 

j-0  J  1  j-0  j  \  L  j-0  J  J, 


(19) 


where  J  represents  the  Jacobian  associated  with  the  pair  of  algebraic  equations.  To  find  the 
stability  boundary  in  question,  one  merely  sets  $  equal  to  zero  in  the  definitions  of  R,  and  1, 
used  in  obtaining  ko  and  k(  in  Equations  (16)  and  (17).  The  result  is  a  boundary  in  the  ko-ki 
parameter  plane  that  is  defined  in  terms  of  system  parameters  dj  [in  the  general  case,  see 
Equation  (2)]  or  a„  bj,  c,  in  the  special  case  of  Equation  (14)  and  the  independent  argument, 
tu,T.  The  independent  argument  is  varied  in  value  between  zero  and  jr  (assuming  that  the 


system  being  designed  possesses  low-pass  filter  characteristics  so  that  only  the“primary  strip” 
need  be  considered),  thereby  defining  the  boundary  in  terms  of  lo,  and  ki.  If  the  coefficients  ajf 
b,.  and  c,  contain  exponential  terms  with  T  appearing  in  the  exponents,  then  each  exponent 
may  be  replaced  by  its  power  series  and  truncated  according  to  the  accuracy  that  is  desired.  If 
the  two  parameters  koand  ki  do  not  appear  linearly  as  expressed  in  Equation  (14),  they  maybe 
solved  for  by  substituting  Equation  (6)  into  the  characteristic  Equation  ( 1 )  and  separating  the 
result  into  the  real  and  imaginary  parts.  This  still  results  in  two  simultaneous  algebraic 
equations  which  may  be  solved  for  the  two  variables,  ko  and  ki,  although  not  as  handily  as 
expressed  in  Equations  (16)  and  (17). 

The  next  two  stability  boundaries  are  those  separating  the  stable  real  roots  from  the 
unstable  ones.  These  comprise  a  mapping  of  the  7=1  and  z=-l  points  from  the  z-plane  onto  the 
selected  parameter  space.  These  are  found  readily  by  substituting  z=l  and  z=-l,  respectively, 
into  the  characteristic  Equation  ( 1 ).  Each  of  the  two  resulting  equations  results  in  a  definition 
of  the  two  real  root  boundaries  in  the  selected  parameter  space. 

The  fourth  boundary  is  a  mapping  of  the  conditions  that  cause  the  two  simultaneous 
algebraic  equations  (used  to  define  the  complex  conjugate  root  boundary)  to  become 
dependent.  In  the  case  where  ko  and  ki  appear  linearly  as  in  Equation  (14),  this  case  is  found  by 
determining  the  conditions  that  cause  the  Jacobian  of  Equation  (18)  to  become  identically 
equal  to  zero. 

For  a  linear  sampled-data  control  system  to  be  stable,  it  is  necessary  that  all  roots  of  the 
characteristic  equation  lie  within  the  unit  circle  on  the  z-plane.  If  the  system  is  low-pass  in 
nature,  the  stable  region  is  bounded  by  the  semi-circle  defined  by  the  upper  half  of  the  unit 
circle  (as  discussed  above)  and  the  singularities  associated  with  z  =  ±1  (as  discussed  above) 
and  J=0  (as  discussed  above).  The  mapping  of  these  boundaries  onto  the  parameter  space  will 
bound  the  stable  region,  if  one  exists.  That  region  is  determined  by  applying  a  shading 
criterion  or  using  a  test  point.1  If  the  Jacobian  is  greater  than  zero,  then  the  stable  region(if  it 
exists)  lies  to  the  left  of  the  complex  conjugate  root  boundary  as  a>„ T  increases;  the  left  side  of 
the  line  is  double  cross-hatched  to  indicate  a  boundary  associated  with  double,  or  complex 
conjugate  roots.  (If  the  Jacobian  is  less  than  zero,  the  stable  region  lies  to  the  right).  Single 
cross-hatching  is  used  on  the  two  contours  associated  with  the  real  roots.  The  side  of  the 
contour  on  which  to  place  the  cross-hatching  is  determined  by  the  requirement  that  cross- 
hatching  be  continuous,  or  on  the  same  side,  of  the  contours  as  the  intersections 
corresponding  to  z=  +i  and  z=  -1  are  approached  along  either  the  complex  root  or  real  root 
stability  boundary.  In  the  unusual  (but  physically  possible)  singular  case  when  J=0,  the 
shading  of  the  resulting  boundary  is  determined  in  a  similar  manner,  i.e.,  the  side  on  which  to 


shade  the  boundary  is  determined  by  the  physical  requirement  that  the  number  of  stable  roots 
must  never  become  less  than  zero  as  one  moves  across  a  stability  boundary.  This  J  =0  case  may 
arise  for  a  particular  frequency,  <o,  or  for  a  particular  combination  of  system  parameters.  The 
latter  situation  is  included  in  the  work  concerning  a  spinning  Skylab.7 

4.  SYSTEM  DYNAMICS 

Once  the  stable  region,  if  it  exists,  has  been  determined  in  the  selected  parameter  space,  the 
dynamic  or  transient  characteristics  of  the  system  can  be  specified  in  terms  of  the  locations  of 
the  roots  of  the  characteristic  equation  (pole  placement).  For  the  complex  conjugate  roots, 
these  locations  are  defined  in  terms  of  damping  ratio  (0  and  system  natural  frequency  (a>n). 
Contours  of  constant  £  are  determined  as  functions  of  cu*  and  T  in  precisely  the  same  manner 
that  the  complex  conjugate  stability  boundary  was  determined  except  that  £  is  not  set  equal  to 
zero. 


The  real  root  locations  corresponding  to  values  of  z  when  6  equals  0°  and  180°  also  may  be 
plotted  on  the  parameter  plane  by  setting  z  equal  to  a  positive  or  negative  real  constant, 
substituting  that  value  into  Equation  (1),  and  solving  for  ki  as  a  function  of  ko.  Each  resulting 
contour  corresponds  to  a  location  of  a  real  root  in  the  z-domain.  When  z=+a 


Z  Y  aj  -  0, 
J-0  i 


and  when  z=  -a. 


(20) 


n/2 
l  Y 
j-0 


2j 


2j 

a 


(n-2)/2 

Z 

j-0 


Y(2j+lj)ct 


(2j+ii) 


0,  n  even 


(n-l)/2 

Z 

j-0 


(Y2ja  Y(2j+l)a 


(2j+l)_ 


0,  n  odd 


(21) 

(22) 


where  a  is  a  positive  real  number.  (The  two  associated  stability  boundaries  may  be  found  by 
setting  a  equal  to  unity).  Values  of  a  corresponding  to  desired  real  root  locations  may  be 
substituted  into  Equations  (20)  through  (22)  to  obtain  parameter  space  contours 
corresponding  to  these  locations. 
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Now  the  dynamic  effect  of  any  chosen  design  point  in  the  parameter  space  may  be  specified 
in  terms  of  £,  <o„.  a.  and  T.  The  analytical  technique  developed  permits  the  designer  to  observe 
the  effect  of  simultaneously  changing  three  control  parameters  and  the  sampling  period.  Most 
existing  conventional  techniques  permit  the  observation  of  the  effect  of  changing  only  one 
control  parameter  and  don't  show  the  effect  of  various  sampling  periods.  If  the  designer  is 
clever,  sometimes  the  system  parameters  may  be  manipulated  within  the  equations  defining 
the  stability  and  dynamic  contours  so  that  more  than  two  parameters  (such  as  ko,  ki)  may  be 
used  to  define  a  parameter  space.  An  example  of  a  three-parameter  space  is  provided  in  the 
sequel. 

Sometimes  it  is  specified  that  the  setting  time  of  the  system  be  less  than  a  prescribed  value. 
This  corresponds  to  requiring  that  the  real  part  of  the  roots  of  the  characteristic  equation  be 
less  than  a  prescribed  negative  real  constant.  A  boundary  corresponding  to  this  requirement 
can  be  drawn  on  the  parameter  plane  by  mapping  a  circle  of  constant  radius  (for  a  chosen 
values  of  o>„,  £,  and  T)  from  the  z-plane  onto  the  parameter  plane.’  Relations  exist  for 
estimating  the  maximum  overshoot  and  peak  time  of  transient  response  when  it  is  valid  to 
assume  a  second  order  system.8  However,  a  simple  estimate  sometimes  can  be  made  merely  by 
looking  at  the  difference  equation  representing  the  system  response,  estimating  when  the 
overshoot  will  occur,  and  plotting  a  corresponding  line  on  the  parameter  plane.  Steady  state 
responses  may  be  found  from  the  open  loop  transfer  function  and  the  assumed  forcing 
functions  in  the  conventional  manner. 

5.  EXAMPLE 

Consider  the  planar  model  that  portrays  controlled  rigid  body  rotational  dynamics.  It  is 
shown  in  block  diagram  form  in  Figure  l.  The  plant  to  be  controlled  is  represented  by  the 
transfer  functions. 

G4(s)  -  4>(s)/Tc(s)  -  1/Jys  (23) 

and 

G5(s)  »  $(s)/$(s)  -  I/s  (24) 


S.  B.C.  Kuo.  Analysis  and  Synthesis  oj  Sampled- Data  Control  Systems.  Prentice-Hall.  New  Jersey.  1963. 


I  I 


Gho(*) 


Gho(») 


Figure  1.  Block  diagram  of  controlled  rigid  body  motion. 


where  s  represents  the  Laplace  operator,  <t>  represents  the  rigid  body  attitude  [in  the  time 
domain,  4(t);  in  the  Laplace  domain,  <£(s)],  Tc  represents  the  commanded  torque,  and  the 
overdot'represents  the  derivative  with  respect  to  time.  The  two  sensors  that  measure  attitude 
( <t> )  and  attitude  rates  (</>)  are  assumed  to  be  perfect  with  unity  transfer  functions.  The  on- 
board  digital  controller  develops  the  commanded  control  torque  from  the  input  states  <t>  and  <p 
and  the  commanded  attitude  signal  Gh0(s)  represents  a  zero-order  hold  in  the  computer 

Gh  (s)  -  <l-e“Ts)/s,  (25) 

o 

where  T  is  the  sample  period  of  the  on-board  digital  computer.  The  PID  control  algorithm  is 
represented  by  Gi(s)  and  Gt(&) 

G^(s)  -  Kp  +  Kj/a  (26) 


where  Kp  and  Ki  are  the  position  and  integral  feedback  gains  and 

V*>  *  *D 


(27) 


where  K»  is  the  derivative  feedback  gain  (more  commonly  called  the  rate  gain).  The  third 
order  characteristic  equation  is  Equation  ( 1 )  with  n=3.  where  the  coefficients  of  z1  are  in  the 
form  of  Equation  (14),  i.e.. 


Y0  *  (a  -  b  +  2c  -  1), 
Yj  -  (“2a  +  2c  +  3), 

Y2  ■  (a  +  b  -  3), 
y3  -  i. 


Modified  gains  a,  b,  c  are  defined  as 


a  =  V/Jv, 


b  5  K  T  /2J  , 
P  V 


c  -  KjT  / 


(28) 

(29) 

(30) 

(31) 


(32) 

(33) 

(34) 


The  stability  boundaries  are  presented  in  terms  of  a.b.c,  and  cu„T  by  observing  that  the  roots 
(in  this  case  three)  of  a  characteristic  equation  representing  a  stable  linear  sampled-data 
control  system  must  all  lie  within  the  unit  circle  in  the  z-plane.  It  is  assumed  that  the  system 
possesses  low-pass  filter  characteristics  so  that  only  the  primary  strip  (corresponding  to  0  <  0 
^  t r)  need  be  considered.  The  stable  region  may  be  defined  by  mapping  its  three  boundaries 
from  the  z-plane  onto  the  a,  b,  c  parameter  space  by  first  considering  parameters  a,  b  to  be  ko, 
k!  in  Equation  (14).  The  boundary  at  z  =  +  l  is  found  by  substituting  that  value  for  z  into  the 
C.E.,  yielding  the  stability  boundary 

c  -  0  (35) 

The  boundary  at  z  =  -l  is  found  by  substituting  that  value  for  z  into  the  C.E.,  yielding 

a  -  2.  (36) 


The  complex  conjugate  root  stability  boundary  may  be  found  by  setting  z  equal  to  cos  a>nT 
+  i  sin  ubT  in  the  C.E.  However,  at  this  point  it  is  convenient  to  use  the  recursive  formulas  of 
Equation  (8)  and  transform  the  C.E.  into  two  algebraic  equations  by  separating  the  real  and 
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imaginary  parts.  Ifthey  are  solved  fora  and  band  the  associated  Jacobian  J,  one  obtains  as  the 
complex  conjugate  root  stability  boundary, 


a  -  Cd  +  B)/(l  -  B)]  c  +  (1  -  B), 
b  -  c  +  (1  -  B), 


J  -  -4  (1  -  B)  <  0  V  0  <  0  <  it.  (39) 

Examination  of  Equation  (39)  reveals  that  the  stability  boundary  associated  with  the  singular 
case  J=0  only  occurs  when  B=l,  which  is  the  already  considered  z=l  case.  Now  the  stability 
boundaries  of  Equations  (35),  (36),  (37)  and  (38)  can  be  plotted  on  a  three-dimensional  plot 
with  axes  a,  b,  c.  A  sketch  of  these  boundaries  is  shown  in  Figure  2.  The  stability  region  is 
found  by  applying  the  “cross-hatching”  rules  or  by  using  one  or  more  test  points  (known  to  be 
stable,  to  lie  on  a  stability  boundary,  or  to  be  unstable).  If  JX)  the  stable  region  (if  it  exists)  lies 
to  the  left  of  the  complex  conjugate  root  stability  boundary  of  Equations  (37)  and  (38),  as  8 
(or  a*»T)  increases;  the  left  side  of  the  line  is  double  cross-hatched  to  indicate  a  boundary 

-a  z=+1  STABILITY 

BOUNDARY  (c=0) 


COMPLEX  CONJUGATE 
ROOT  BOUNDARY 

{a  -  b  -  2c  [1  -  (b-c)"1] } 


2=  -1  STABILITY 
BOUNDARY  (a»2) 


Figure  2.  Stability  boundaries  in  3-dimensions. 


associated  with  double,  or  complex  conjugate,  roots.  If  J<0.  as  in  this  case,  the  stable  region 
lies  to  the  right  and  the  right  side  of  the  line  is  double  cross-hatched.  Single  cross-hatching  is 
used  on  the  contours  of  Equations  (35)  and  (36)  associated  with  single  root  boundaries.  The 
side  of  the  contour  on  which  to  place  the  cross-hatching  is  determined  by  the  requirement  that 
cross-hatching  be  continuous,  or  on  the  same  side,  of  the  contours  as  the  intersections 
corresponding  to  z=+l  (0=0)  and  z=-l  (0=jt)  are  approached  along  either  the  complex  root 
or  real  root  stability  boundary.  Using  this  criterion,  the  stable  region  is  found  to  be  bounded 
by  the  c  =  0  plane,  the  a  =  2  plane,  and  the  curve  specified  by  Equations(37)  and  (38).  The  latter 
may  be  solved  for  a  as  a  function  of  b.  resulting  in  the  simple  expression 

a*b  -2c  Cl-(b-c)  *  D  .  (40) 

The  limiting  values  fora  and  b  may  be  found  by  letting  d  approach  0  and  n,  i.e.,z=±l.  For  the 
case  of  d  —  0,  represent  B  by  the  series  expansion 

B  -  cos  0  -  l-(02/2!)+(04/4!)-(06/6!)+...  (41) 

Then  using  Equations  (37),  (38).  and  (41),  one  obtains 


lim  b  -  c 

0-0  (42) 

lim  a  .  -  lim  c  [(4/02)-(2/3)+(02/36)-. . . .]  -  « 

0-0  0-0  (43) 

For  c=0.  Equation  <42)  still  holds  (recognize  ki=0  is  really  the  case  one  is  interested  in).  Then, 
using  Equations  (34)  and  (37),  one  obtains 


lim  a  -  lim  (k_03/4J  ui3)  [2-(02/2!)+(04/4!) - ] 

0-0  0-0,  V  , 

(0V2 1 )- (0*/4  !)+••• 

-  lim  (kj0/2J  w3)[2-(202/6)+(04/72) - ]-0 


(44) 


For  the  limiting  case  of  9  —  jr: 


lin  a  -  lim  {C(l+B)c/(l-B)3+(j-B)}  -  2 
©  •+  IT  ©  -*■  IT  (45) 


and 


lim  b  ■  lim  (c+l-B)  ■  c+2 
0  IT  0  ■*  IT 


(46) 


The  results  of  the  foregoing  are  sketched  graphically  in  two  dimensions  on  Figure  3,  for  the 
cases  where  c  =  0  and  when  c  0. 

The  design  technique  employs  the  pole  placement  approach.  It  is  based  on  the  premise  that 
the  dynamic  behavior  of  the  system  is  related  closely  to  the  location  of  the  roots  of  its 
associated  characteristic  equation.  The  method  shows  both  analytically  and  graphically  the 
direct  correlation  between  these  roots  and  the  control  gains  and  sample  period  of  the 
controller.  The  design  technique  then  involves  the  specification  of  the  characteristic  equation 
root  locations  and  the  subsequent  determination  of  the  control  system  gains  and  sample 
period  needed  to  attain  these  locations.  The  control  system  designer  then  must  determine  the 
system  response  resulting  from  using  these  numerical  values  and  assess  its  adequacy.  If  it  is  not 
adequate,  he  usually  relies  on  his  experience  to  relocate  the  roots  to  improve  the  response  in 
the  manner  desired  for  his  particular  system  (i.e.,  faster  settling  time,  lower  peak  overshoot, 
etc.).  It  is  assumed  that  one  wishes  the  pair  of  complex  conjugate  poles  of  the  C.E.  to  dominate 
the  dynamic  response  of  the  system.  This  response  will  then  be  modified  by  locating  the  third 
(real)  root. 

First  consider  the  pair  of  complex  conjugate  roots.  One  may  use  Equations  (16)  -  (19), 
substituting  a,  b  for  ko,  k>,  respectively  [also  using  Equations  (28)  through  (31)].  Upon 
examining  the  nature  of  the  equations  for  a,  b,  one  sees  that  they  are  functions  of  r,  B,  and  c,  or 
(finally)  functions  of  £,  9,  c.  For  each  value  of  (,  b  versus  a  may  be  determined  analytically  and 
plotted  graphically  as  a  function  of  9  for  a  given  value  of  c.  This  may  be  repeated  for  a  number 
of  values  of  c.  Now  the  value  of  9  may  be  selected  (for  selected  values  of  c)  that  best  meets 
whatever  design  criteria  one  chooses.  In  effect,  selecting  9  and  c  will  also  choose  the  real  root 
location.  This  may  be  shown  by  designating  the  real  root  of  the  C.E.  as  zr,  where 

zR  =  6  -  e  R  (47) 
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Figure  3.  Stability  boundaries  In  2-dimensions. 


!  :  I  i. 


and  or  represents  the  damping  factor.  One  may  substitute  Equation  (47)  into  C.E.  (I)  and 
solve  for  b,  obtaining 

b  -  [(l-5)/l+6)]a+[2/(l-6)  c-  ( 1-6) 2/ (1+6)]  (48) 

One  sees  that  Equation  (48)  yields  a  straight  line  contour  of  b  vs.  a  for  a  given  value  of  c  and  a 
given  value  of  6.  Thus,  for  a  given  value  of  c.  where  the  6-contour  of  Equation  (48)  crosses  the 
(•contour  of  Equations  (37)  and  (38),  one  obtains  a  value  of  6  for  the  given  value  of  c  and  those 
specified  values  of  (  and  6.  From  Equation  ( I ),  c  may  be  determined  as  a  function  of  the 
characteristic  equation  root  locations: 


c  -  (a2  +  $2  -  2a  +  1) (1  -  6)/4, 


where  a  and  0  are  the  real  and  imaginary  parts,  respectively,  of  the  pair  of  complex  conjugate 
roots  of  characteristic  Equation  (1).  For  a  specified  location  (a,0)  of  the  pair  of  complex 
conjugate  roots,  a  value  for  c  is  established  by  setting  a  value  for  6  (location  of  the  real  root). 

An  additional  constraint  is  imposed  by  Shannon’s  Sampling  Theorem: 

«  <  «8/2,  (50) 

where  at  is  defined  in  Equation  (S)  and  uhI2it  is  the  sampling  frequency. 


ui  -  2tt/T. 

s 


From  Equations  (SO)  and  (51),  this  constraint  may  be  restated  as 


0  <  ir  . 


The  stable  region  portrayed  in  Figures  2  and  3  may  be  shown  for  several  selected  values  of  c 
on  Figure  4  using  Equations  (16),  (17)  and  (48).  As  c  increases  in  value,  the  stable  region 
shrinks  in  size  and  the  {  contours  lie  further  to  the  right.  When  c=0.6,  the  (-  1  ,VTcontour  lies 
entirely  to  the  right  (and  outside)  of  the  stable  region.  Finally  when  c^2,  the  stable  region 
disappears  entirely  (for  a  geometric  explanation,  see  Figure  2). 


4 


The  application  of  the  foregoing  technique  will  now  be  summarized. 

Based  on  desired  system  response  characteristics,  tentatively  select  locations  for  the  three 
roots  of  the  characteristic  equation.  Although  the  zeros  of  the  system  characteristic  equation 
also  affect  the  dynamics,  placement  of  the  poles  will  dominantly  affect  the  dynamics. 

•  A  root  location  criteria  might  be  to  select  numerical  values  for  £  and  a*,. 

•  Determine  computer  design  constraints.  One  might  be  to  select  as  large  a  value  of  T  as 
possible  to  minimize  on-board  memory  capacity  requirements. 

4  Then  qualitatively  determine  how  much  integral  gain  is  desired  for  the  types  of 
disturbance  inputs  expected.  This  gives  an  indication  of  the  value  of  c  to  use. 

Using  the  value  of  c  tentatively  established  in  the  above  paragraph,  plot  the  stability 
boundaries  in  the  b  vs  a  parameter  plane,  using  Equations  (36)  through  (38).  Then  plot  the 
desired  {-contour  (as  a  function  of  argument  tonT)  on  the  b  vs.  a  parameter  plane,  using 
Equations  ( 16)  and  (17). 

Using  various  values  of  6.  plot  6-contours  on  the  same  plot,  using  Equation  (48). 

Find  an  intersection  of  a  6-contour  which  gives  a’ desired  value  for  to„T  (and  hence  T,  since 
a>n  has  already  been  specified)  and  6.  In  general,  if  6  —  0.  its  effects  on  the  system  dynamics  is 
small  compared  to  the  effect  of  the  pair  of  complex  conjugate  poles  (placed  by  o>BT,  £). 

Check  the  response  of  the  system  using  the  values  of  a.  b,  c.  and  a>nT  (and  hence  K.p,  KLi.  Kd, 
and  T)  associated  with  the  selected  intersection.  If  it  is  unsatisfactory,  reiterate  the  above 
procedure,  selecting  another  value  of  c. 

If  the  design  procedure  were  merely  to  specify  the  three  root  locations  (such  as  in  terms  of  a, 
/3,  6,  a>„T.  6),  iteration  would  not  be  required.  Unfortunately,  in  design  practice  one  must 
exercise  engineering  judgement  in  selecting  the  three  root  locations,  giving  rise  to  iterative 
design  procedures  such  as  the  one  outlined  above.  However,  another  possibility  does  exist.  If 
the  real  root  is  placed  near  the  origin  of  the  z-plane  (6  —  0),  and  if  the  pair  of  complex 
conjugate  roots  are  located  near  the  unit  circle  in  the  z-plane  (i.e..  specifying  £,  and  wnT).  the 
latter  two  roots  will  dominate  the  system  response.  Then  the  already  developed  tools  for 
specifying  the  system  response  (such  as  determining  explicitly  the  maximum  overshoot  and 
the  settling  time)  of  second-order  systems  may  be  applied  handily.  The  foregoing  procedure 


may  be  amplified  through  the  application  of  realistic  numerical  values  for  the  system 
parameters. 

A'S'.ime  ?  digital  on-board  controller,  and  further  assume  that  it  is  desired  to  have  a 
controller  natural  frequency  (o>n)  of  six  rad  s,  and  that  the  moment  of  inertia  (J.)  about  the 
single  axis  is  4.3  x  104  kg.m:.  Assume  that  integral  control  is  desired  to  drive  to  zero  the  effect 
of  constant  input  disturbances.  This  means  a  non-zero  value  of  c  is  desired.  If  a  value  of  {  of 
1  jyjl  is  selected,  and  if  the  value  of  6  is  desired  to  be  kept  as  close  to  zero  as  possible,  it  is 
implied  from  Figure  4c  that,  for  5  >0,  no  intersections  of  the  {-contour  occur  when  c>  0.2. 
Hence,  a  value  of  c  =  0.1  is  chosen.  The  smallest  value  of  6  whose  corresponding  6-line 
intersects  the  {-contour,  with  a  reasonable  factor  of  safety,  is  6  =  0.4.  The  value  of  WT  =  1 .2  is 
on  the  {-contour  at  the  intersection,  yielding  a  value  of  0.20s  for  T.  Corresponding  values  of  a 
and  b  are  1.41  and  0.68,  respectively.  This  gives  a  set  of  gain  values  of  1.14  .x  I07  N-m/rad,  1.67 
x  10"  N  m  rad-s,  and  2.35  x  10ft  N  m  s/  rad  for  K.,,.  Ki.  and  Ki>.  respectively.  If  response  studies 
show  that  a  larger  value  of  integral  gain  is  deemed  necessary,  a  larger  value  of  c  (and 
consequently  a  smaller  value  of  {)  would  be  selected. 

6.  CONCLUSIONS 

An  analytical  method  for  portraying  stability  regions  in  a  selected  parameter  space  has  been 
shown  for  a  digital  system.  The  method  requires  that  the  system  characteristic  equation  be 
available  and  expressed  in  the  complex  z-domain.  It  also  is  possible  to  apply  pole  placement  to 
ootain  desired  dynamic  characteristics  using  this  modified  parameter  space  technique.  The 
advantage  of  the  technique  over  existing  classical  sampled-data  methods  is  that  the  stability 
and  dynamic  response  characteristics  are  expressed  in  terms  of  several  (rather  than  merely 
one)  selected  parameters.  Also,  the  sampling  period.  T.  need  not  be  expressed  numerically 
before  the  design  technique  begins,  giving  the  system  designer  one  more  degree  of  freedom. 
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It  is  possible  to  obtain  algebraic  recursive  relations  for  the  real  and  imaginary  parts  of  the 
complex  variable.  /.  when  it  is  raised  to  the  jlh  powerlj  is  a  possible  integer).  From  Equation 
(3),  z  is  defined  as  the  vector  (or  complex  variable), 

z  ■  eTs  ■  r  r  cos  0  +t  r  sin  0.  (Al) 

Using  Equation  (3)  in  Equation  (Ai),  i.e„  letting  B  =  cos  0,  leads  to 

z  -  rB  +  i  r  vi-B2.  (A2) 

The  real  and  imaginary  parts  of  z1  may  be  defined  as  R,  and  1„  respectively,  i.e., 

zJ  -  R  +  i  I  ,  (A3) 

If  j=l,  Equation  (A3)  becomes  identical  with  Equations  (AI)  and  (A2),  leading  to 

-  rB  (A4R) 

and 

II  -  ry/l-B*  .  (A4I) 

If  j=0,  then  Equation  ( A3)  becomes 

z°  -  1  -  Rq  +  i  IQ  ,  (A5) 

or 

Rq  -  1  (A6R) 

and 

I0  -  0.  (A6I) 
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if  one  squares  the  left  and  right  hand  sides  of  Equation  (Al  i,  one  may  obtain 


z2  .  r2  ei29»  r2  (  cos  29  +  i  sin  29) 

*  r2  (  2  cos29  -  1)  +  12  r2  sin  9  cos  9 

2  2 
■  2  r  cos  9  (cos  9  +  i  sin  9)  -  r 

*  2rBz  -  r2.  (A7) 

One  may  now  multiply  each  term  of  Equation  (A7)  by  z**1  to  yield 


lt+1  „  _  k  2  k-1 
z  »  2rBz  -  r  z 


(AS) 


Substitution  of  Equation  (A3)  into  Equation  (A8),  letting  j  *  k-1,  k,  and  k+I,  leads  to 


Vi  *l\n-  2rt«ic  +  1  V  -  r  <\-i+  Irk-1> 


(A9) 


Separately  equating  the  real  and  imaginary  parts,  respectively,  of  Equation  (A9)  leads  to  the 
two  equations, 

Rk+1  -  2rBR^  -  r2Rk_1  (A10R) 

and 

Vi  -  2riS  -  '2  Vi  • 

One  may  observe  the  identical  forms  of  Equations  (A  10)  and  write  them  both  as  a  single 
equation, 

Vi  ■  2rBjtk  -  'Vi  •  <ai» 


where  Xk  may  be  used  to  represent  either  Rk  or  lk.  Of  course  the  dummy  index  k  may  be 
exchanged  withj,  yielding  Equation  (8)  of  the  text.  If  one  knows  the  two  values  X*  and  Xk-i  of 
Equation  (All),  one  may  determine  the  third  value  Xk*i.  Since  two  sets  of  values  of  Xk  are 


known,  i.e..  for  k  -  0  and  k  =  i  [see  Equations  (A4)  and  (A6)],  as  many  values  of  Xk  as  are 
needed  for  the  particular  problem  at  hand  can  be  determined  recursively  from  Equation 
(All). 
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THIS  PROGRAM  IS  STORED  UNDER  MI  OPR 

THE  PROGRAM  CALCULATES  THE  VALUES  OF  PI  AND  P2  IN  THE 
PARAMETER  PLANE  METHOD.  THE  STABILITY  BOUBARIES  APE 
CALCULATED  IF  ZETA  =  0.  IN  ADDITION  THE  PROGRAM  WILL 
CALCULATE  THE  VALUES  OF  PI  AND  P£  FOR  a  LINEAR  SET  OF 

equations. 


OPTION  BASE  e 
DEFAULT  ON 
FLOAT  4 

D I M  C  (  2  0  >  ,  S  <  2  O  >  ,  A  <  2  0  )  ,  E  <'  2  0  >  ,  F  <  2  0  >  ,  DU  <20  1 ,  M  m  U0,3>,E  %  4  <.2  0  )  ,  T  h  *  <  2  6  > 
D I M  Tent l  <  20  > ,  T *i»2 < 20  > ,  T 2a >  ,  T*m < 20 >  ,  Rr  < 20 >  ,  Ir< 20 ) ,  R<  20  > ,  D 1  < 20 > 
DIM  Zb < 20)  ,  Wnb< 20  > ,  Ph i  < 20  "> ,  A  i  pha< 20 > 


THE  USER  IS  GIVEN  THE  OPTION  OF  HAVING  DATA  PLOTTED  AND/OR 
HAVING  THE  RESULTS  PRINTED.  THE  PRINTER  IS  DESIGNATED  BY  EITHER 
"7,6"  FOR  THE  THERMAL  PRINTER  OR  '‘16"  FOR  THE  CRT" 


DISP  “  This  is  r-hs  MI  CPF  program.  PI  ease  prsss  CONTINUE" 

PAUSE 

L  INPUT  "IS  DATA  TO  BE  PLOTTED?  Y  OR  N  ?".Qp» 

Ltyp«-1  !  SETS  COUNTER  FOR  PLOTTER 

IF  <  QpS*  " Y“ )  OR  02p**“y")  THEN  GOSUB  2790 

L INPUT  "DO  YOU  WISH  A  HARD  COPY  OF  THE  PROGRAM  RESULTS?  Y/N“ fT*st* 
PRINTER  IS  7,6 

IF  <T*st#»“N")  OR  <T*st  #«“«“.>  THEN  PRINTER  IS  16 
PRINT  LINOS) 

L  I HPUT  " S , 2 ,  OR  W  TRANSFORM? " , T 1  * 

<  THE  USER  HERE  HAS  THE  OPTION  OF  SELECTING  EITHER:  2  ,  W  OR  3  TRANSFORM 
! 

IF  T I S  *  "  S "  THEN  PRINT  "S  TRANSFORM  METHOD" 

IF  T 1 S  =  "  2  "  THEN  PRINT  "2  TRANSFORM  METHOD" 

IF  TII='Z"  THEN  530 
IF  T I  *  * " S "  THEN  540 
BEEP 

DISP  "  ONLY  2  OR  $  ALLOWED,  "  ;  T  1  * ;  "  ,  WAS  INPUT,  PRESS  CONTINUE  TO  GO 

PAUSE 
GOTO  400 


INPUT  “WHAT  IS  THE  ORDER  OF  THE  CHARACTERISTIC  EQUATION?". N 

PRINT  USING  56a; N 

IMAGE  "THE  SYSTEM  ORDER  IS  *" , DD 

PRINT  LIN'.  1> 


530  1  , 

530  ! 

$00  ! 

•510  ' 

620  ! 

620  REDIM  C<N>  ,S<N>,  ft<N>  ,  B  N  >  ,  F<N>  ,  Dd<N>  ,  Rr  ■  N 

*40  INPUT  "INPUT  ELEMENTS  OF  A  MATRIX,  fU0x,A‘ l 

$50  PRINT  "ELEMENTS  OF  A  MATRIX  A <  0>  ,  A <  1  > ,  A 2  :> ,  .  .  .  "  ,  A * >  ; 

660  IF  <01 p*=“Y" )  OR  <Glp*="y"  >  THEN  600 

670  INPUT  "INPUT  ELEMENTS  OF  E  MATRIX,  B< 0 <, B< 1 B< * > 

©SO  PRINT  "ELEMENTS  OF  B  MATRIX  B<  0>  ,  B<  1>  ,  B< 2  B<  *> ; 

6S0  IF  <  0 1  p#*“  Y  “  >  OR  <  0 1  p  *  *  “  y  “  >  THEN  710 

700  INPUT  "INPUT  ELEMENTS  OF  F  MATRIX.  F  O  - ,  F  <  l  j  F  <  * 

710  PRINT  "ELEMENTS  OF  F  MATRIX  F<  OJ  ,  F  < 1  >.  F  <2>,...n,  F<  4  )  • 

7  30  ! 

740  ' 

7 SO  t  THE  value  of  zeta  -DAMPING  RATIO..'  is  input  here  for  use  in  the 
770  !  PARAMETER  SPACE  METHOD. 

7  SO  ! 

7'?0  ! 

300  INPUT  "FOR  THE  PARAMETER  PLANE  CALCULATIONS,  WHAT  IS  VALUE  OF  2ETA^".Z«t* 
310  IF  Zet  a> 1  THEN  GOTO  14S0 
3 2 O  PRINT  LIN < 4 > 

330  PRINT  "ZETA  Zeta 

340  Ltyp=Ltyp*l 

350  Paii=0  !  SET  VALUES 

3S0  IF  Z  *  t  a  >  1  THEN  GOTO  14S0 

370  !  IF  ZETA  IS  EQUAL  TO  ZERO  THE 

330  IF  <T1#«“Z">  AND  <Z«t**0)  THEN  G03UB  2$20  !  PROGRAM  CALCULATES  THE 
330  •  STABILITY  BOUNDARIES  FOR  Z-  +  -1 

3O0  ! 

310  !  j 

320  !  INPUT  THE  VALUES  OF  INITIAL  AND  FINAL  OMEGA  TO  BE  USED.  IN  ADDITION 

330  !  IF  THE  NUMBER  OF  STEPS  TO  BE  USED  IS  EQUAL  TO  ZERO,  THEN  THE  PROGRAM  I 

340  i  ASSUMES  THAT  THE  USER  WISHES  TO  SKIP  THE  SECTION  ON  THE  PARAMETER-  \ 

350  '  SPACE  AND  SKIPS  DIRECTLY  TO  THE  OPTION  TO  EITHER  CALCULATE  THE  ROOTS  ] 

3$0  !  j 

370  i 

330  INPUT  "INPUT  INITIAL  OMEGA,  MAXIMUM  OMEGA,  NO  OF  STEPS" , W  Hit , Maxw, Nwnt  j 

330  PRINT  "INITIAL  OMEGA  »" ? Mi nt ; TAB < 35 > ; " F INAL  OMEGA  *";M*xw  j 

1000  IMAGE  ,  "NUMBER  OF  STEPS*" , MDDD, 

1010  PRINT  USING  1O00; Hunt  j 

1020  IF  W  i n t * O  THEN  Wint  =  lE-25  j 

1030  IF  Nwnt  <  *0  THEN  14$0  j 

1O40  ! 

1O$0  ! 

1070  !  THE  CALCULATIONS  ASSOCIATED  WITH  THE  PARAMETER -SAPCE  DESIGN  FOLLOW 
1030  1  ! 

IO30  ! 

1100  1  ' 


1110  IF  C  T  I  #  =  "  3  "  3  OR  <  T  1  *  =  "  j  "  >  THEN  Th*r  i»P  I-ASN<$QRt  1 -2*t  v>-2>  > 

1120  FOR  I*0  TO  Nwm 

1  1  30  Wnt  *W  i  nt  + 1  *  (Maxw*W  i  n»-  >  Hwnt 

1140  GOSUB  2250 

1150  GOSUB  2110 

1150  IF  X*Q  THEN  GOTO  1200 

1170  IF  T  1$  *  "  2  "  THEN  PRINT  "  UriT  =  “  i  Unt  ;  THE  >.  20  >  ;  "  BET  =  '* ;  X;  TAB  <  40  •  "p  l  =  •• ;  pc.  0  >  ;  TAB'  5 
> ? “ P2* “ ; P  < t  > 

1 100  IF  T  1  #  =  "  S  “  THEN  PRINT  “Hn*"  :  Mm  :  TAB<20>  J  “D£T»"  j  JCj  TAB<40>  .*  "Pl»“j  P.  0  >  j  TftBfe-0 
5 “P2»";P<1> 

1100  GOTO  12S0 

1200  IF  T 1 #* " S "  THEN  PRINT  “Wn  *"  :  Writ 
1210  IF  T]#*"2"  THEN  PRINT  “WnT  *";Wr,t 

1220  PRINT  "THE  DETERMINANT  OF  THE  COEFFICIENT  MATRIX  IS  EQUAL  TO  ZERO" 

1230  PRINT  "MATRIX  AB  =“  i  IF  SO  THEN  SKIP  THE  NORMAL  ROUTINE 

1240  MAT  PRINT  Ab 
1250  PRINT 

1250  PRINT  "MATRIX  Fp  *" 

1270  MAT  PRINT  Fp 

1200  IF  ■:  Qpi  *  “  V "  )  OR  CQpJ="y">  THEN  GOSUB  2820 
1200  NEXT  I 

1300  ! 

1320  ! 

1330  !  A  SERIES  OF  CONTROL  STATEMENTS  FOLLOW  THAT  GIVE  THE  USER  VARIOUS  OPTIONS 
1340  !  TO  CHANGE  CERTAIN  VAR  ISLES. 

1350  ! 

1350  ! 

1370  IF  •;  Qpf*  "  Y  "  )  OR  (  Qp#*“y  “  :>  THEN  GOSUB  2338 

1.380  L  INPUT  "CHANGE  OMEGA  RANGE-,  2ETA?,  OR  F  MATRIX?  THEN  INPUT  W,Z,F  OR  N0",Q1* 
1390  IF  <Q1**"N0")  OR  03 1  f  *  “  N  “  3  OR  <Qli  =  "riO"  >  OR  <Ql*  =  "ri">  THEN  1458 
1400  IF  Ql#»“W"  THEN  9S0 
1410  IF  01#*  "2"  THEN  800 

1420  IF  (Q1#*"F">  AND  (01 p#*"Y" >  THEN  1380 
1430  IF  Ql#*"F“  THEN  700 

1440  IF  ( Q 1 1 <  >  " W "  )  OR  (01  #s  >  " Z "  >  OR  ( 0 1  #<  >  “ F "  .)  THEN  BEEP 
1450  IF  (Qt#<>“ W" >  OR  (01 f< >"2" >  OR  03 1 $< > “F" >  THEN  1380 
1450  LINPUT  "DO  WISH  TO  CALCULATE  THE  CONTOURS  FOR  REAL  ROOTS?", Rr# 

1470  IF  R r J * " Y "  THEN  GOTO  3180 
1430  • 

14  90  ! 

1500  '  THE  USER  IS  GIVEN  THE  OPTION  OF  CALCULATING  THE  COEFFICIENTS  OF  THE 
1510  '  CHARACTERISTIC  EQUATION. 

1520  ! 

1530  1 

1540  LINPUT  "CALCULATE  COEFFICIENTS  OF  THE  CHAR.  EQU.  OR  ROOTS?  Y  N?" , T* 

1550  IF  ( T # * ” H "  OR  <T**“n“  >  THEN  1958 

1558  INPUT  "INPUT  PI , P2" , Cpl , Cp2  !  PI  IS  ASSOCIATED  WITH  A  AND  P2  WITH  £ 

1570  MAT  T*m 1 *A* (Cp  1  ) 

1580  MAT  T*m2*8*(Cp2> 

1590  MAT  T*m3*T«» 1+T*m2 
1500  MAT  T*m*T*m3+F 

1510  PRINT  "COEFFICIENTS  OF  THE  CHARACTERISTIC  EQUATION  WITH" 

1520  PRINT  "Pl*"j  Cpl ;  TAB',  30> ;  "P2*";  Cp2 
1530  PRINT 


1640  PRINT  "C<0>  C<l>*S  ♦  C<2>*S-'2  +  C<3>*3'3  ♦  ..." 

1650  FOR  1=0  TO  N 

1660  PRINT  USING  1  670  ;  I ,  T*i»  ■  I' 

1 670  IMAGE  “C< “ , DD, “ >=“ . MZ. 6DE 
1 6  8  0  T  •=  m  1  <  I  ■  =  0 
1 630  NEXT  I 

1700  DISP  "WAIT  FOR  SILJAK" 

1710  1 

1720  !  Oil)  for  SILJAK 

1730  ! 

1  740  CALL  Si  I  j*k  <  N,  ?*#><*.>  ,Tfr»l^',.  ©000001 ,  .  0O00001 , 50.  Rr  ■  *  • .  I  •  *>'.  • 

1750  IF  T  I  J  =  "  S  “  THEN  PRINT  L  I N  (  2  i  ,  “ROOTS :  “  ,  L  I N  •,  1  >  ,  SPA  •;  8  ,  "  REAL  “  ,  "  IMAGINAR 
DAMPING  RATIO","  Wrn”  ,  L IN < 2 > 

1760  IF  T  1  *  *  "  Z  "  THEN  PRINT  LIN  <  2  >  ,  "ROOTS  :  "  ,  L  I N  ■.  1  >  ,  SPA  ■.  3  .• ,  “  REAL  "  ,  "  IMAGINAR 

DAMPING  RATIO","  WrlT",LIN<2> 

1770  FOR  1=1  TO  N 

1780  IF  ABS  < I r  < I >  H  *  1 £- 1 1  THEN  Ir<I>=0 

1730  R  t  * m  =  <  Rr  <.  I  > -'2+ 1 r  <  I  >  A2 > .  5 

1  SOO  S  t  ■=  m=Rr  <  I  >  "Rt.  fin 

1310  IF  I r < I >  =0  THEN  St*m=l 

1320  IF  T  1  f  =  "  S  "  THEN  GOTO  1330 

1330  Sgmt*LOG<Rt 

1340  Omgt  =ACS  <  St-  *(ti ) 

1 350  Om*g*n t  *SQR  <  Sgmt  -2+0 r,,. 3 1  ^2  > 

1360  Z*  t  ant  *-Sgmt  'Omeg  ant 
1370  GOTO  1300* 

1330  Om#g*nt »Rt «« 

1330  Z*t in t»-St* m 

1  SOO  PR  I  NT  US  I  NG  1 330 ;  Rr  <  I  >  ,  I  r  <  I  >  ,  2s  ♦.  ant ,  Omgirn 

1310  NEXT  I 

1320  PRINT  L I N < 5 > 

1 330  IMAGE  ?X , MZ. 6DE , 5X, MZ. 6DE , 6X, MZ . 6DE , SX , MZ . 6DE 
1340  GOTO  1380 

1350  IF  ( 0  p  f  *  "  Y "  >  OR  <Qp$*Hy“*)  THEN  GOSUB  3000 
I960  GOTO  200 
1370  1 
1330  ! 

1330  ! 

2000  ! 

2010  1 
2020  ! 

2030  ! 

2040  ! 

205O  1 
206O  ! 

2070  ! 

2038  !  _ 

2030  !  CALCULATE  A  ,B  AND  F  MATRICES  AT  EACH  OMEGA 
2100  ! 

2110  Ab 0 , 0 )  =DOT <  A ,  C > 

2120  Ab  <  0 , l ) *DOT<  B , C  > 

2130  Ab<  1 , 0)  =  DOT<.  A,  S) 

2140  AtKl, l)»DOT<B,S) 

2150  Fp(0)«-DOT(F,C) 

2160  F  p  < 1 > » - D  0 T  < F , S  > 

2170  DIM  D <  1  ,  1  >  ,  At.  •,  1 ,  1  >  ,  P <  1  > ,  F p <  1  > 

2130  MAT  D= I NV 1 Ab  > 

213 0  MAT  P  =  D * F p 

2200  :i  =  DET 

2210  RETURN 
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25© 

FOR  K»0  TO  N 

2*0 

IF  <T1*-“Z">  OR  < 

270 

Rk=Wnt ■ K 

2  SO 

Rrgk  *Th* t 

230 

GOTO  2320 

300 

Rk«EXP<-2«*  -a* Writ* 

310 

Rr-gk  =  K.*Wnt  #«.  1  -Z*». 

32© 

C<K>*Rk*COS<. Argk  > 

330 

S'. K  >  =  R k* S  I N Argk  :> 

340 

NEXT  K 

350 

RETURN 

360 

1 

370 

f 

330 

1 

33© 

| - - - - - 

400 

| 

4  10 

( 

■  T  !  f  =  "  -  "  THEN  230© 


THIS  SUBPROGRAM  WILL  SET  THE 
VALUES  OF  THE  C  AND  S  MATRICES 


FOR  K  =  1  TO  N, 
ZETA  AND  WN*T 


GIVEN  VALUES  OF 


THE  SUBROUTINE  TO  CALCULATE  CONSTANT  ZETA  CONTOURS 


430 

i 

440 

i 

450 

FOR  1=0  TO  N 

460 

Dd< 1 >  =  Dk AI 

470 

NEXT  I 

430 

E  =  DOT  < A, Dd> 

430 

L=DOT<B,Dd> 

500 

M=BOT  < F , Dd> 

510 

IF  LOO  THEN  PRINT 

“  P2=  “ :  -E  •  L ;  “  *P  1  +  " 

;  -CML 

520 

IF  <L=0>  AND  ' E< >0 • 

THEN  PRINT  "Pi  ="; 

-M,  E 

53© 

IF  <L«©>  AND  <E«0J 

THEN  PRINT  "COEPFIC 

I  ENTS 

540 

PRINT 

550 

PRINT  "THE  AXES  INTERCEPTS  ARE  AS  FOLLOW 

S :  " 

560 

IF  LOO  THEN  PRINT  “AT  P1  =  0  ,  P2  =  "  i  -M-’L 

57© 

IF  EO0  THEN  PRINT  "AND  AT  P2*0  .  Pl  =  “{- 

M-  E 

530 

IF  L=0  THEN  PRINT  "PI 

=  CONSTANT  ="j-M.'E 

530 

IF  E=0  THEN  PRINT  "P2 

=  CONSTANT  =" j-M/L 

60O 

PRINT  LINC2) 

*  1  © 

RETURN 

620 

<  l  > 

63S 

PRINT  "THE  STABILITY  BOUNDARIES  FOP 

Z  =  +  1 

Dk  =- 1 

640 

PRINT  "FOR  Z=-l“ 

650 

GOSUB  Elm 

660 

PRINT  L  I N  <  2  '< 

670 

Dk  =  1 

630 

PRINT  "FOR  2  =  + 1 " 

630 

GOSUB  Elm 

7O0 

PRINT  L I N  <  2  > 

71© 

RETURN 

EQUAL 


•1  APE  AS  FOLLOWS “ , L 


i  j  j  I  !  j  M  i  i  m  1  i  t  j  i  i  t  i  i  I  i  i  i  i  i  i  i  i  m  i  i  i  i  f  fi  f  f  f  |  !  !  !  !  !  M  !  1  1  !  !  f  1  f  ?  M 
THE  FOLLOWING  SECTIONS  BETWEEN  THE  "MM!”  ARE  FOR  PLOTTING  1 


2790 
2800 
2810 
2S20 
2S30 
2S40 
2050 
2S00 
2370 
2830 
2800 
2000 
2010 
2020 
2030 
2040 
205@ 
2000 
2070 
2080 
2000 
3000 
3010 
3020 
3030 
3  0  4  0 
3050 
3000 
3070 
3080 
3000 
3100 
3110 
3120 
3130 
3140 
3150 
3100 
3170 
3130 
3  1 00 
3200 
3210 
3220 
3230 
3240 


ASSIGN  #1  TO  “ F'DAT  A" 

READ  #1,4 
RETURN 

IF  R r S <  > “  Y “  THEN  PRINT  #  1 ;  P c 0 >  ,  P <  1  > ,  0 ,  Set  a,  END 
IF  Rr*="Y"  THEN  PRINT  #  1 ;  P < 0  ■ ,  P ■:  1  • .  0 ,  1 +Dk  ,  END 
IF  Rr*="Y"  THEN  GOTO  2020 
IF  Pass*0  THEN  Mi  ;<*Max=PC0> 

IF  Pass*©  THEN  Miy*M*y*P( 1 > 

P  as  s  s  1 

IF  P < 0 > <  =M i x  THEN  Mix»PC0> 

IF  P<0>  >«Max  THEN  Max=Pc0> 

IF  P(lX*Hiy  THEN  Miy=P<l> 

IF  PC  1 >  >  =  May  THEN  H*y«P(l) 

RETURN 


PR  I  NT  tt 1 ; P  <  0  > , P  Cl > , 1 , Lt  yp , END 
REDIM  MmCLtyp, 3> 

M  M  <  L  t  y  p ,  0  >  *  M  i  x 


MmCLtyp,  1  '.i  =Max 
MmcLtyp, 2>*Mi y 
MmCLtyp,  3)=t1iy 
RETURN 

PRINTER  IS  10 

PRINT  “The  program  is  now  over  and  t  he  data  is  on  the  tape. 
PRINT 

PRINT  "If  you  would  like  to  plot  the  data" 

PRINT  “then  you  must  toad  the  plot  routine  entitled  PLOTER. 
IF  CTest*="y“  >  OR  C  Test  *  =  “  Y "  >  THEN  PRINTER  IS  7,0 

REDIM  Mm<Ltyp,3> 

PRINT  #1 ,  1  ; Ltyp, Mm<# > 

RETURN 


II 


i* 


;  i  i  t  i  I  !  ;  i  i  |  !  m  i  |  <  i  i  |  ;  m  i  !  i  i  i  i  i  i  i  i  i  i  i  i  i  i  i  i  !  i  i  i  i  i  i  i  i  i  i  i  i  i  i  i  i  |  i  !  i  ;i  ; 


THIS  PART  OF  THE  PROGRAM  CALCULATES  THE  CONSTANT  SETA  CONTOURS 


!  CALCULATE  the  boundaries  at  SETA  >1.0 

PRINT  "FOR  REAL  ROOTS  OF  THE  CHARACTERISTIC  EQUATION” , LINCl > 
INPUT  "WHAT  IS  THE  VALUE  OF  DELTA?”, Dk 
PRINT  "DELTA  *”;Dk 
GOSUB  Elm 
PRINT 

IF  < QpJ  =  ” V " >  OR  C  Qpi= " y " >  THEN  GOSUB  3320 


3250 

3260 

3270 

3230 

3200 

3300 

3310 

3320 

3330 

3340 

3350 

3360 

3370 


GOTO  1460 


Subroutine  to  prepare  values  for  plotting  of  the  constant 
seta  contours  for  real  roots. 


IF  LB0  THEN  GOTO  3420 

INPUT  "INPUT  MINIMUM  s,  MAXIMUM  VALUES  OF  P  1  ”  ,  P 1  m  i  ,  P 1  ma 
P<0>»Plmi 

PC  1  >«-<E*P<0>  +M>.'L 
GOSUB  2320 
PC0>»Plma 
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3380  P<1 >«-<E*P<0>+M>/L 

3390  GOSUB  2820 

3400  GOSUB  2930 

3410  GOTO  3498 

3420  INPUT  "THE  CONTOUR  IS  VERT  I CHL , THUS  INPUT  MAX, MIN  VALUES  OF  P 

3430  P(0)*-M,'E 

3440  P  <  1  ) *G 

3450  GOSUB  2820 

3480  P< 1 ) *H 

3470  GOSUB  2820 

3430  GOSUB  2930 

3490  RETURN 

3500  I 

3510  ! 

3520  ! 

3530  END 
3540  ! 

3550  1 

3580  ! 

3570  ! 

3580  ! 

3590  ! 

3800  1 

3810  ! 

3820  ! 

3830  SUB  Si  1  j  4k  <  N ,  Rc  o*f  <  * ) ,  Icoof <*> .  T  o  l  a,  To  1  f ,  1 1.  m  ix ,  Rroot  <  *'• ,  I  root  <  *  >  > 
3840  BIN  Xa i 1 j ik < Os N> , V* i 1 j*k < 0 8 NJ 
3850  ! 

3880  ' 

3870  ! 

3880  ! 

3890  !  #»*  POLYNOMIAL  ROOTF IHDER  USING  SILJAK'S  METHOD. 

3700  !  4«****4**»*»*#»************»»*#'#»****«*»4'******* 

3  710  ! 

3720  ! 

3730  !  ***  BAD  DATA  CHECK. 

3740  Biddt **<N< »0)  OR  <Tola<*0>  OR  <To1t'<*Q>  OR  <ltm*x<«0> 

3750  IF  Baddta=0  THEN  3890 

3780  !  ***  PRINT  ERROR  MESSAGE  AND  PAUSE. 

3770  !  ***  USER  MAY  CORRECT  DATA  AND  CONTINUE. 

3780  PRINT  LIN<2>, "ERROR  IN  SUBPROGRAM  $ i 1 j ak . " 

3790  PRINT  “ N»  *  j  N, " To  1 ** “ j  To  1  * 

3800  PRINT  "Tol  f  ■  "  ;  T  o  1  ,  "Umix*"  ;  I t »*x, L INC 2> 

3810  PAUSE 
3820  GOTO  3740 
3330  ! 

3840  ! 

3850  !  ***  BEGIN  SUBPROGRAM. 

3880  !  ***  INITIALIZE  LOCAL  VARIABLES. 

3370  !  ***  RrootC*)  AND  Iroot<*>  HOLD  THE  ROOT  APPROXIMATIONS. 

3330  •  ***  INITIALIZE  Rroot <»>  AND  Iroot<*>  TO  9.999999E99. 

3390  MAT  Rroot » < 9 . 999999E99 > 

3900  MAT  I  root  =  <  9. 999999E99> 

3910  ! 

3920  ! 

3930  Nn*N 

3940  !  ***  SPECIAL  CHECK  FOR  N-l. 

3950  IF  N* 1  THEN  5080 
3980  ! 

3970  ! 
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398©  !  #**  Xs  i  1  j  ak  <  *  >  AND  Y*i  1.jak<*>  HOLD  THE  SILJAK  COEFFICIENTS. 

3990  Y*Yi i  1 j ak  < 1 > *Xa i t  j  ak • 8  >  *  1 
4O80  X»Xa j I jak  < 1 1 
4810  Yj i 1 jak <0>»L*0 
4820  ' 

4038  ! 

4048  '  ***  BRANCH  TO  COMPUTATION  OF  SILJAK  COEFFICIENTS. 

4058  GOSUB  Si  1  j  ak 

4868  !  •***  G  GETS  FORMER  VALUE  OF  F,  INITIALIZE  AND  INCREMENT 
4878  !  *•**  ITERATION  COUNTER. 

4038  G*F 
4090  M  *  Q  *  P  ■  8 
4108  L*L* 1 

4110  !  #4*  Z»  <  DU'  DX  >  '2*  C  DV  DX  >  A2 . 

4138  FOR  K-l  TO  N 

4130  P»P*K#<Rco#f <K>*Xt i l jak<k-l>- Ieo*f <K> #Y*i 1 j ak <K-t > > 

4140  Q»Q*K#  <  Rco*f<K>*Y*1  1  jak  CK  - 1  >  +  Icoef  <  K  )  #X*  i  I  j*k  <K- 1  )  > 

4150  NEXT  K 
41o8  Z*P*P+Q*Q 

4  170  !  ###  De  l  t.  ax  AND  D«U*y  ARE  THE  RESPECTIVE  CHANGES  IN  X  AND  CHANGES  IN 
4138  D«1 tax*-CU*P#V»Q>  '2 
4198  D*  1 1  ay®  CU*Q-V#P  >  -'2 

4208  !  ***  INCREMENT  SUCCESSIVE  QUARTERING  COUNTER. 

4218  M«M+ 1 
4220  ! 

4238  ! 

4240  •  #»#  NEW  ROOT  APPROXIMATIONS  X, Y  LOADED  INTO  Xiiljak<l>  AND  Y#f 1  J4k< 1> 
4258  Xi  i  1  .j ak  C  l  >*X  +  D*  1 1  ax 
4268  Y* i 1 j  ak  C 1 ) ■ Y  +  De It  ay 

4270  (  ***  RECOMPUTE  SILJAK  COEFFICIENTS. 

4238  GOSUB  Si  1 j  ak 
4298  ! 

4380  ! 

4310  !  ***  IF  NEW  ERROR  ESTIMATE  GREATER  THAN  OLD,  QUARTER  THE  SIZE 
4320  !  ***  OF  Del  tax  AND  D* 1 1 ay . 

4338  IF  F  >*G  THEN  4538 

4348  !  #*#  CHECK  IF  INCREMENTS  ARE  SMALL  ENOUGH  TO  SATISFY  THE 
4350  !  ***  STOPPING  CONDITIONS. 

4368  IF  CABS c  De  1 1  ax  .X  To  1  a>  AND  < ABSv  D*  1 1  ay  XTo  1  a>  THEN  4790 

4378  !  ###  CHECK  IF  MAXIMUM  NUMBER  OF  ITERATIONS  HAS  BEEN  EXCEEDED. 

4330  <  ***  IF  SO,  PRINT  ERROR  MESSAGE  AND  PAUSE. 

4390  IF  L > 1 1  max  THEN  4788 
4488  ! 

4410  ! 

4420  !  ***  ITERATE  AGAIN. 

4438  X  »  X  s  i 1 jak  C 1 ) 

4440  Y»Y»itjakCl> 

4450  GOTO  4080 
4468  ! 

4470  ! 

4438  <  ##«  QUARTER  THE  INTERVAL  SIZE  AND  TRY  AGAIN. 

4490  <  ***  MAXIMUM  #  OF  QUARTER  I NGS  HAS  BEEN  SET  AT  20. 

4508  !  ***  IF  THIS  IS  EXCEEDED,  THEN  A  FINAL  TEST  IS  MADE  ON  U  AND  V  TO 

4518  •  *#*  SEE  IF  THE  FUNCTIONAL  TOLERANCE  IS  SATISFIED.  IF  NOT,  AN 

4528  *  *##  ERROR  MESSAGE  IS  PRINTED  AND  THE  PROGRAM  EXITS. 

4538  IF  M>28  THEN  4570 
4540  D* 1 t ax*D* 1 t axx4 
4558  D*  1 1  av*D#  1 1  ay^'4 
4560  GOTO  4210 
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4570  IF  <ABS<U:><=To]  f  >  AND  <ABS<V>  OTol  t\>  THEN  479© 

458©  ! 

4590  ! 

480O  !  ***  PRINT  ERROR  MESSAGE  AND  PAUSE. 

4810  PRINT  LIN<2>  ,  "ERROR  IN  SUBPROGRAM  Sill  it.  .  " 

4820  PRINT  "THE  INTERVAL  SIZE  HAS  BEEN  QUARTERED  2©  TIMES  AND  " 

4830  PRINT  "THE  TOLERANCE  FOR  FUNCTIONAL  EVALUATIONS  IS  STILL  HOT  MET. 
4840  PRINT  "Tol  t'«"  ;  To!  f,  “U»“;U,  "  V  =  "  ;  V  ,  L  I N  <  2  > 

4850  PAUSE 
4880  ! 

4873  ! 

4880  !  ***  MAXIMUM  #  OF  ITERATIONS  HAS  BEEN  EXCEEDED. 

4890  ■  ***  PRINT  ERROR  MESSAGE  AND  PAUSE. 

4700  PRINT  LIN<2>, “ERROR  IN  SUBROUTINE  Siljak." 

4710  PRINT  "MAXIMUM  «  OF  ITERATIONS  HAS  BEEN  EXCEEDED." 

4720  PRINT  "La " ; L ,  "It  max*"  ;  1 1  max , L I N( 2 '> 

4730  PAUSE 
4740  GOTO  4390 
4750  ! 

4780  ! 

4770  !  4»*  ROOT  FOUND.  STORE  COMPUTED  ROOT  IN  Rroot <*>  AND  Iroott*) 
4780  !  ***  ARRAY  ELEMENT. 

4790  Rroot  <N>*Xsi  l  jak  <  1  > 

4800  Iroot < N > *Ys i 1 j ak < 1 > 

4810  ! 

4820  ! 

4830  !  ***  INTIALIZE  VARIABLES  FOR  SYNTHETIC  DIVISION. 

4840  A*Rcoef<N) 

4850  BaIco*f<N> 

4880  Rco*t'<N)*Ico*t'<H>*0 
4870  X*Xsi 1  jak  <  1  > 

4380  Y*Y*  i  1  jak  <  1  > 

4890  !  ***  SYNTHETIC  DIVISION  TO  CALCULATE  NEW  COEr  FI  CENTS  Rco*t'<*> 
4900  !  ***  AND  Ico«f<*>. 

4910  FOR  K-N-l  TO  0  STEP  -t 
4920  CaRco«f<K> 

4930  DaIco*t'<K) 

4940  U^RcosC^K+l ) 

4950  V*Ico*f <K+1 > 

4980  Rc  0*f < K  >  *A*X*U-Y*V 

4970  Icoaf  OOaB+X*V  +  Y*u 

4980  A*C 

4990  B*D 

5000  NEXT  K 

5010  N»N- 1 

5020  !  ***  REDUCE  NUMBER  OF  COEFFICIENTS  AND  BEGIN  AGAIN. 

5030  IF  NOl  THEN  3990 

5040  !  ***  SINCE  DEGREE  OF  RESULTANT  POLYNOMIAL  IS  ONE, 

5050  !  ***  COMPUTE  FINAL  ROOT  ALGEBRAICALLY. 

5060  A*Rco*t'<0> 

5070  U*Rco*F<U 
5080  B»lcoef<0> 

5090  V* Ico*f  < 1 > 

5100  T »U*U+V*V 

5110  Rroot  < 1 >»-<A*U+B*V>xT 

5120  Iroot<l>*<A»V-U*B>/T 

5130  N*Nn 

5140  SUBEXIT 
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5150  ! 

5130  ! 

5170  !  ***  SUBROUTINE  TO  COMPUTE  SILJfiK  COEFFICIENTS. 

5130  S I ! jak s  2»Xs i 1  jak  < 1 )*X* i 1 j  ak  < 1 >+Y*i 1 jak  <l>*Y*i 1 jak < 1 
5130  T*2*X» i 1 j ak < 1 ) 

5200  FOR  K»0  TO  N-2 

5210  Xai 1 jak  <K*2)«T»X*  i  1  jak  CK+l  > -Z*Xs  i  1  jak  <IO 
5220  Ya 1 1 jak  <K+2> =T*Ys i 1 jak  <K+ 1 > -2*Ys i 1 jak  <K  > 

5230  NEXT  K 
5240  U*V“0 
5250  FOR  K*0  TO  N 

5230  U*U+Rcoa*'<K>*X»i  1  jak  <K>-Ico**'<K>*Y*i  t  jak 

5270  V*V*Reo*f<IO»Ys1  1  jak  <K>*Ico*f  <K>*Xai  1  jakOO 

5230  NEXT  K 

5230  F»IJ*U+V*V 

5300  RETURN 

5310  SUBEND 


*  4  * 


10 

20 

30 

40 

70 


1  00 
t  10 
120 
1  30 
140 

150 
i  SO 
170 
ISO 
1  SO 

r  y  y 
210 

240 

250 


THIS  PROGRAM  CALCULATES  THE  RESPONSE  USING  THE 
CROSS  MULTIPLICATION  TECHNIQUE 


IT  IS  STORED  UNDER  "CROSS" 


rt**Tir**- rff»*--t*r4*»4*****if^rff**t**f* 


,  M  fh  ■  .3 ,  10: 


REM 
REM  * 

REM  ■* 

REM  4 
REM  * 

REM  4 
REM  4 
REM  4 

REM  4  *  *  •*•*.*  4  * 

OPTION  BASE  0 

DIM  B  v  25  > ,  A  •  25  >,R<  200  >  ,  Cc  c  25  • ,  Rr  '.25  > ,  C  '  200  . 

L  =  -l 

L  INPUT  "DO  YOU  WISH  A  HARD  COPY  OF  DATA  AND  RESULTS’" , One  * 
PRINTER  IS  IS 

IF  Qhc**"Y"  THEN  PRINTER  IS  7,6 

L INPUT  "DO  YOU  WISH  TO  HAVE  THE  RESULTS  PLOTTED  ■  V  N  ;■  "  ,  Q p I 
IF  Qp*»"N"  THEN  GOTO  200 
ASSIGN  #1  TO  "PDATA" 

READ  #1,3 

INPUT  "INPUT  THE  ORDER  OF  ThE  NUMERATOR,  DENOM I NATOR " , N , D 
RED  IN  A ' N  > , £ i D  > 

PRINT  USING  270; N,D 

IMAGE  "NUMBER  OF  STEPS  FOR  RESPONSE  *" . ODD 
INPUT  "INPUT  SAMPLE  INTERVAL", T 
PRINT  "SAMPLE  INTERVAL  =";T 


260 

GOSUB  Input. 

270 

IMAGE  "ORDER  OF 

NUMERATOR  IS  “ , DD, ION, "ORDER  OF  DENOMINATOR  IS  " , D D 

2  SO 

GOTO  540 

2 S  0 

Inp:  i  Subroutine 

t  o  c  al  cu  1  at  &  i  nput  f  u n c  t  i  on 

3  0  0 

LI NPUT  "DO 

YOU  WISH  UNIT  STEFCU) , SINCS) , RAMPCR > , OR  USER  DEFINED'  D  ■ 

N  AS 

INPUT? " , Qr* 

3  1  0 

IF  Qr t  =  " IJ “ 

THEN 

GOSUB  Step 

320 

IF  0 r # ■ " S " 

THEN 

GOSUB  Sin 

3  30 

IF  Qri»"R" 

THEN 

GOSUB  Ramp 

34  0 

IF  Qr4*"D" 

THEN 

GOSUB  U**r 

350 

IF  <Qr*<>" 

U“>  AND  <  Qr#<  >  "  3  '*  >  A  N  D  <  Q  r  *  <  >  "  R  “  ,•  A  N  D  C  0  r  t  <  >  "  D "  >  THEN  G  0  T  0 

360 

GOTO  410 

370 

BEEP 

330 

DISP 

"ONLY 

U,  S,  R,  OR  D  ALLOWED" 

3  SO 

WAIT 

3500 

4  0  0 

GOTO 

300 

4  1  0 

RETURN 

4  20 

SUBEND 

4  30 

! 

440 

1 

450 

1 

4  0 

I nput. : PRINT 

L  I N  C  2 

>, “CALCULATIONS  MADE  ACCORDING  TO  FOLLOWING  EQUATION 

TION 

II 

4  70 

PRINT  LINC 

1  >  ,  " 

C  <  2  )  •■■■  R  c  2  )  «  N U M E R A T 0 R  D E N 0 M  I  NAT 0 R  " 

4S0 

PRINT  “ 

NUMERATOR  =  AC0)#2  '0  +  A<  J  >*2*'l  ♦  AC2)#2  '2  .  .  .  ACI 

4*0 

PRINT  " 

DENOMINATOR  =  BC0)*2A0  4  B<1  >*2M  +  BC2)#2'S2  .  .  .  BC 

500 

RETURN 

510 

SUBEND 

520 
5  30 
540 
550 
560 
570 
5S0 

5  SO 
600 
610 
620 

6  4  y 


D  I SP  "  TO  N»"  ?  N ;  "  I NPUT  NUMERATOR  COEFF I C I  ENTS  A  < 0  ;■ ,  A  ■  1  > ,  A  <  2  >  ,  .  .  .  A  ■  N  >  "  ; 

INPUT  A  <  *  > 

PRINT  L  I N  <  '3  :< .  "NUMERATOR  COEFFICIENTS  A  CO  > ,  AO  >  ,  A<  2  >  ,  .  .  .  A<N.>  "  ,  AC#)  ; 

L I NPUT  " ANY  CHANGES’  c Y • N ■ " . Qc h if 

IF  Qcha**"N"  THEN  GOTO  630 

INPUT  "COEFFICIENT  NUMBER  I , A C I  • " , I , A •  I  ■ 

PRINT  USING  610; I, AC  I  > 

IMAGE  "AC",k,">  «  "  ,  K. 

GOTO  570 

DISP  "TO  D  =  “  :  D  s  "  INPUT  DEMCMImAtOR  COEFFICIENTS  B  0  .  E  ■  1  • ,  B  ■  2  ■ ,  .  .  .  £  •  D  • 
input  b-4. 


650 
660 
670 
660 
690 
700 
710 
720 
7  3  0 
740 
750 
760 
770 
760 
760 
800 
610 
620 
690 
640 
650 
660 
870 
660 
880 
900 
910 
920 
9  30 
940 
950 
960 
970 
980 
990 
1000 
1010 
1020 
1030 
1040 
1050 
1060 
1070 
1050 
1  096 
1  1  00 
1110 
1120 
1  1  30 
1  140 
1150 
1  1 60 
1170 
1  180 
1  1  90 
1200 
1210 
1220 
12  30 
1240 
1250 
1260 
1270 
1280 
1290 
1  300 


PRINT  "DENOMINATOR  COEFFICIENTS  Bt0>  ,  B<  1 .  B<  2)  .  .  .  .  B<  N>  "  ,  B<* '> » 

L INPUT  “ANY  CHANGES?  < Y/N> ", Och4* 

IP  0  c  h  * “ N "  THEN  GOTO  720 

INPUT  "COEFFICIENT  NUMBER  I ,  B ■  I  -  "  .  I ,  B -  I  > 

PRINT  USING  610;  I,  BU> 

IMAGE  " B < " , K ,  “  >  *  “ , K 
GOTO  660 
Bd-B- D> 

RED IM  B- D-l  >  ,  Rr • N  > . Cc •  D-  1  • 

L  *L  *  1 

INPUT  “INPUT  NUMBER  OF  STEPS  FOP  RESPONSE  ",  Ns t e p 
PRINT  USING  2  30:  Ns t  ep 
REDIM  R  Nst  ep  > 

GOSUB  I  rip 

FOR  K«0  TO  Nstep 

FOR  J«0  TO  N 

IF  t>J-D<0  THEN  Rr  •  .1  >  *0 

IF  K*J-D< 0  THEN  GOTO  840 

Rr  >.  J  >  *R  f  K* J-D  > 

next  ; 

FOR  J*0  TO  D-l 

IF  K+.T-D  <*0  THEN  Cc  ’  J  f  *0 

IF  K*J-D<»0  THEN  GOTO  690' 

Cc • J>«C<K*J-D> 

NEXT  J 

C‘,K>*I  Bd*<  DOT-  A.Pr /'-DOT-  B,  Cc  •> 

IF  K  »  O  THEN  -1 1  V  a  M  iy  =  C  0  > 

IF  Mi y > C < K >  THEN  Mi y * C ( K > 

I  F  M*y<  C  '•  K  >  THEN  Miy*C  <  K  ■' 

IF  Qp**"Y"  THEN  GOSUB  Tips 
NEXT  K 

IF  Qp*»"N"  THEN  GOTO  1020 

PRINT  #1  :  Nstep*T,  RCNi  t  «p  >  ,  O'  Nstep  .  1  .  L+l  ,  END 

Mm  <  O ,  L;-*0 

M  m  < 1 , L>»Nstep*T 

Mmt  2, L)*Mi  y 

Mm<  3, L  > *  M  iy 

GOSUB  Output 

L  INPUT  "CHANGE  THE  INPUT  DATA-.:  IV'  NUMBER  OF  SAMPLES ■: N  >  *  OR  3T0P<.S>  ?"  ,  Qi  p* 

IF  uipi»',N"  THEN  GOTO  750 

IF  Qi p**" I “  THEN  GOTO  20O 

IF  Q i p t < > " 8 ”  THEN  BEEP 

IF  Q  i  p  f  <  >  "  $ "  THEN  GOTO  10  30 

IF  Qp J* “ THEN  GOSUB  Final 

END 

Step:  !  Routine  to  jet  "ilues  For  P  when  input  is  unit  step 
PRINT  "THE  INPUT  IS  A  UNIT  STEP  FUNCTION" 

MAT  R-CON 

RETURN 

SUBEND 


Sin:!  Routine  to  set  wilues  to  R  when  input  is  SIN  function 
PRINT  "INPUT  FUNCTION  IS  SIN  WAVE" 

INPUT  "INPUT  FREQUENCY  IN  PAD- SEC , AMPL 1  TUBE  ", W , X 
PRINT  "FREQUENCY  -PAD  SEC  '  IS  *";W 
PRINT  "AMPLITUDE  OF  SINE  WAVE  = " : X 
FOR  K * 0  TO  N s  t  e  p 
P-  K.  ;■  ■  X * S  I N  -  W*K*T  • 

next  k: 

RETURN 

SUBEND 


G-2 


*•  *  riWj-l, 


I  ' 

l  V 


11 

I- 1 


1  3 1  0 

R  amp :  1  R out  i  f> 

to  Sit  "%lu*J  of  P 

h  *  ft 

1  300 

PRINT 

"INPUT  IS  RAMP" 

1  3  SO 

INPUT 

“IN P  U T  3  L  0  P  E  0  F  PAM P 

■  PEP 

1  340 

PRINT 

"SLOPE  IS  *"iW 

1350 

FOR  K* 

0 

TO  Nit*p 

1360 

ROO«W 

*T 

1  370 

NEXT  K 

1330 

RETURN 

1  300 

SUBEND 

1400 

i 

14  10 

i 

1400 

t 

1  430 

U  ■■  *  r  s  1  R  u  t  t  r.  *  t  o  s  *  * 

of  R  for  us* 

r  d« 

1440 

PRINT  "THE  INPUT 

IS  A  USER  DEFINED 

FUNC 

14  50 

FOR  K * 0  TO  Hit 

*p 

14  6  O 

FIXED  0 

1470 

DISP  "INPUT  R' 

5 

K ;  "  >  "  ; 

1430 

INPUT  R'.K> 

14  00 

NEXT  K 

1500 

PRINT  "INPUT  FOR 

USER  DEFINED  FUNCTION 

1510 

STANDARD 

1500 

PRINT  R  < * > 

15  30 

RETURN 

1540 

SUBEND 

1550 

0  u  t  p  u  t  i  1  Rout  t  fi*  to 

pr 

j  ft  t  t  hi  d  it  i  f  *  S  U  1  t 

£ 

c' 1  r  t  i  t  i  p;« 


I  560  PRINT  L I N <.3 '■ 

1570  PRINT  ’*  K"  ?  SPA<SJ ;  "TIME"  jSPflai)  :"  INPUT"  j  SPAk  1 
150O  FOR  K*0  TO  Nst*p 

1500  PRINT  USING  l«®0jK,K*T,ft<K>,C<K> 

1600  IMAGE  ,  DDD5X  ,  MD  .  DDBE ,  lOX.MD.DDBE,  10X.  HD.  BODE 
1610  NEXT  K 
1600  RETURN 
1630  SUBEND 

1640  Tip* :  !  PLACE  THE  DATA  ON  TAPE  FOR  PLOTTING 

1650  PRINT  #1}K#T,R<K>,C<K>,0,L*1,END 

1660  RETURN 

1670  SUBEND 

1600  1 

1600 

1700  ' 

1710  Ftr-.il  FINAL  DATA  ON  TAPE  ROUTINE 
1700  PRINT.  #1  .  1  ?L 
1730  FOR  J*0  TO  L 

174  0  PRINT  #1  ;  Mm<  0.  J  • ,  Mm  ■  1  ,  I  ,  M»<  0.  J  :• ,  Mm>.  3,  Jh  .  END 
1750  NEXT  .T 
1760  RETURN 
1770  SUBEND 


'RESPONSE' 


CONTROL  DYNAMICS  COMPANY 

701  COfVLETT  DRIVE  -  SUITE  2 
HUNTSVILLE,  ALABAMA  3SS02 

April  2,  1980 

Mr.  Janies  McLean 
US  Army  Missile  Command 
Redstone  Arsenal ,  AL  35809 

ATTN:  DRSMI-RGN 
Dear  Jim: 

In  accordance  with  your  instructions  I  participated  in  the  PERSHING  II 
Autopilot  Design  Review  at  The  Martin  Company  in  Orlando  on  31  March  1980. 

At  the  conclusion  of  the  review,  I  made  a  number  of  comments  on  my  observations. 

The  Martin  Company  representatives  requested  that  Mr.  Jesse  Armstrong  (PERSHING 
Project  Office)  make  all  such  comments  available  in  writing.  Mr.  Armstrong 
agreed.  My  comments  are  summarized  in  this  letter.  I  have  not  included  the 
comments  made  by  Mr.  Ed.  Herbert. 

1.  General . 

The  autopilot  work  that  was  reported  upon  was  excellent,  as  were  the 
presentations.  The  participants  appeared  extremely  capable.  The  results  of 
the  linear  analyses  conducted  appeared  quite  satisfactory,  although  examples 
of  associated  analyses  were  not  presented.  The  review  was  presented  at  a  very 
fast  pace.  It  is  recommended  that  the  next  such  review  be  conducted  at  a  pace 
more  conducive  to  technical  interplay  between  the  presenters  and  the  reviewers. 

This  would  be  enhanced  by  providing  all  reviewers  with  copies  of  the  material  to 
be  reviewed  --  perhaps  a  week  before  the  review. 

2.  Nonlinearities. 

My  recollection  is  that  only  intentionally  included  nonlinearities  (e.g. 
the  saturation-type  nonlinearity  following  each  digital  "integrator")  were  re¬ 
cognized  as  existing.  It  was  assumed  that  the  various  system  states  would  lie 
within  these  saturation  bounds,  hence  precluding  the  application  of  a  linear 
analysis.  I  am  concerned  with  the  possibility  of  other  nonlinearities  existing 
within  the  controlled  missile  plant,  actuators,  and/or  sensors  and  their  dynamic 
effect.  I  do  not  believe  this  possibility  has  been  given  much  attention.  The 
analysis  of  nonlinear  effects  on  a  digitally-controlled  plant  is  difficult,  time- 
consuming,  and  tedious.  Still,  it  has  been  used  in  the  oast  to  uncover  unexplained 


TELEPHONE 
205  /  MI-1342 


simulation  dynamics  and  proven  its  worth.  On  other  programs  I  have  always 
found  rewarding,  if  difficult. 


3.  Quantization. 

Quantization  is  a  special  nonlinearity  associated  with  digital  con¬ 
trollers.  It's  dynamic  and  accuracy  effects  apparently  have  not  been  investi¬ 
gated. 

4.  Modeling. 

During  the  presentation,  it  was  brought  to  light  that  the  dynamics  of 
all  elements  of  the  system  might  not  be  known  (e.g.  engine  nozzle  character¬ 
istics).  One  technique  to  include  their  effects  would  be  to  incorporate  the 
actual  hardware  in  the  loop. 

5-  Tolerance  analysis. 

In  the  planned  tolerance  analysis,  significant  attention  should  be 
devoted  to  the  variation  of  system  parameters.  Two  areas  in  particular  should 
receive  attention: 

(1)  Bending  characteristics.  During  the  presentation,  the  statement 
was  made  that  the  NASTRAN  program  is  producing  bending  parameters  with  an 
accuracy  of  5%.  My  experience  with  NASTRAN  leads  me  to  question  this  statement. 
The  significance  of  accurate  bending  data  becomes  particularly  crucial  when  one 
attempts  to  use  "notch  filters"  to  compensate  for  body  bending  dynamics. 

(2)  Pole-zero  cancellation.  Several  times  during  the  presentation, 
pole-zero  cancellation  techniques  were  alluded  to.  Pole-zero  cancellation  is 
very  sensitive  to  not  knowing  orecise  numerical  values  of  system  parameters. 

I  would  be  wary  of  such  an  approach. 

6.  Computer  capacity. 

Ouring  the  presentation  it  was  stated  that  the  computer  is  already 
filled  to  70%  capacity.  Typically,  at  this  stage  of  development,  computer 
requirements  would  be  expected  to  grow  more  than  30%.  Consequently,  efforts 
should  be  implemented  to  keep  this  growth  to  a  minimum  and  perhaps  to  invest¬ 
igate  the  possibility  of  decreasing  (through  innovation)  computer  requirements 
now. 

As  a  result  of  the  foregoing  observations,  the  following  recommendations 
are  offered  for  consideration. 

(I)  Provide  copies  of  the  proposed  presentation  material  to  all  reviewers 
approximately  one  week  before  the  next  technical  review. 


1 


(2)  Investigate  the  possible  existence  of  dynamically  significant 
nonlinearities.  Predict  their  dynamic  effect  by  mathematical  analysis  to 

be  verified  by  computer  simulation.  Where  this  is  deemed  impractical,  include 
in  the  hybrid  simulation  the  actual  hardware  that  gives  rise  to  the  nonlinearity. 

(3)  Investigate  the  sensitivity  of  the  control  system  to  variations 

in  numerical  values  of  bending  parameters.  View  the  5%  accuracy  figure  associated 
with  the  NASTRAN  program  with  skepticism.  This  should  have  a  strong  effect 
on  the  planned  "notch  filters." 

(4)  Investigate  analytically  the  dynamic  (and  accuracy)  effects  of 
digital  computer  quantization. 

(5)  Investigate  the  sensitivity  of  the  control  system  to  realistic 
variations  in  numerical  values  of  system  parameters  involved  in  any  pole-zero 
cancellation  techniques. 

(6)  Initiate  action  to  decrease  computer  capacity  requirements  and/or 
strongly  limit  future  additional  demands  upon  the  onboard  computer. 


I  would  appreciate  your  reviewing  this  letter  and  forwarding  it  to  Mr.  Jesse 
Arsmtrong.  Thank  you. 


cc  --  Dr.  H.  Pastrick,  DRSMI-RGN 
Mr.  E.  Herbert,  DRSMI-RGN 
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